Métodos iterativos para problemas inversos

Transcript

1. Determinação de raízes de equações Prof. Paulo R. G. Bordoni UFRJ Métodos iterativos para problemas inversos Prof. Paulo R. G. Bordoni
2. Na aula anterior estudamos métodos diretos para resolver problemas inversos. O nome se justifica porque eles fornecem a solução através de um número finito de operações, para construir diretamente a solução, usualmente através de fórmulas.
3. Na aula de hoje os Mestres iniciarão os métodos iterativos para a resolução problemas inversos, em particular de equações
4. Indo direto ao ponto: Assuma, por enquanto, que equações são expressões do tipo = 0. Sim, e que desejamos obter valores de para os quais () = 0.
5. Sejam mais precisos. Quero saber: 1. Qual é o domínio da função , isto é, em que conjunto vou procurar ? 2. E em qual conjunto toma valores ? Seremos dúbios, Loirinha: Algumas vezes nos referiremos a como função e outras como expressão () da equação.
6. Sua pergunta faz sentido, Loirinha. Uma ideia é começar com funções de domínio real e a valores reais: : ⊆ ℝ → ℝ, ↦ = (). Nesse caso e seriam números reais e os valores de que tornam = 0 seriam raízes da equação.
7. Mas, ampliando um pouco o horizonte, poderíamos pensar em funções : ⊆ ℝ2 → ℝ, (, ) ↦ = (, ). Nesse caso o gráfico de (, ) é uma superfície em ℝ3 e a equação , = 0 é a interseção dessa superfície com o plano . Mais precisamente, a solução da equação , = 0 (se existir) é uma curva no plano .
8. Equações do tipo , = representam a fértil ideia de curvas de nível (o nível ), amplamente utilizadas em geociências. Se a é linear, isto é, , = + então a superfície , = é um plano. Portanto + = 0 é a reta onde esse plano corta o plano .
9. Ampliando ainda mais o horizonte, poderíamos pensar em funções : ⊆ ℝ2 → ℝ2, , ↦ 1 2 = , = 1 (, ) 2 (, ) . Nesse caso o gráfico de (, ) está em ℝ2 × ℝ2 uma hiper superfície em ℝ4. Entretanto, o gráfico da 1ª função coordenada 1 = 1 (, ), define uma superfície em ℝ3. Idem para a 2ª função coordenada. 2 = 2 (, ) define outra superfície em ℝ3.
10. Elipse Parábola Hipérbole Circunferência 2 + 2 = 1 2 = 4 2 2 + 2 2 = 1 2 2 − 2 2 = 1 Assim, a solução de um sistema de 2 equações a 2 incógnitas (caso exista) é uma curva em ℝ3. Na “Introdução à MatPlotLib” vimos o caso em que uma superfície é um cone e a outra um plano.
11. No caso de duas equações lineares ቊ 1 , = 11 + 12 2 , = 21 + 22 o sistema fica ቊ 11 + 12 = 1 21 + 22 = 2 Um sistema linear de 2 equações a 2 incógnitas cuja solução, se existir, é a interseção de dois planos: a reta em que eles se cortam.
12. Voltemos agora ao caso mais simples: Obter a solução de uma equação = 0.
13. Uma hipótese imprescindível para existência de solução é que a função : , → ℝ seja contínua em , . Vamos assumir também que a equação já está na forma padrão = 0. Senão podemos reescrevê-la como = − ().
14. Uma classe muito ampla de algoritmos para resolver essas equações são os métodos iterativos. Iterar é sinônimo de repetir. Os mais conhecidos são o da bisseção, o de Newton-Raphson, o do ponto-fixo e algumas variantes deles.
15. Na realidade, nos métodos iterativos, vamos achar aproximações r* para a raiz r. Em outras palavras acharemos números ∗ tais que − ∗ < para ε bem pequeno. A bem da verdade épsilons minúsculos, por exemplo = 10−7, 10−10, 10−15, etc.
16. Nesse caso, 2 − 1.4142135623730951 < 0.5 × 10−16. É o aspecto prático, Mestra, é difícil precisarmos conhecer 2 com mais que umas 5 ou 6 casas decimais após a vírgula. Contei 16 casas na resposta de Python: 2 = 1.4142135623730951 …
17. A aproximação r* para a raiz r será obtida a partir da sequência , , , ⋯ , , ⋯ quando ela converge para r. Quando ela converge, para algum k suficientemente grande, teremos − < . Então escolhemos ∗ = .
18. Mestres, imploro por algo mais prático! Surfista, a parte conceitual tem que preceder a prática. Deixe de ser apressado!
19. Ok Surfista, vou explicar o método da bisseção, também conhecido como método do meio-intervalo.
20. Assumiremos que : [, ] → ℝ é uma função calculável através da expressão () e que é contínua em [a, b]. é o valor de f em . ()
21. () () () Pelo Teorema do Valor Intermediário, se f possui sinais contrários nos extremos a e b então existe pelo menos um ponto (, ) tal que = 0. Dito de outra forma, a equação = 0 possui pelo menos uma raiz em (a, b). Num ponto , () poderá ser positiva, negativa ou nula:
22. = = Como 1º passo do processo, definimos o extremo esquerdo fazendo ↤ (leia: receber a) e o extremo direito fazendo ↤ . Descoberto um intervalo (, ) no qual f possui sinais contrários, podemos começar o processo iterativo. • Em programação: ↤ • Em matemática: =
23. O 2º passo é calcular o ponto médio do intervalo, ↤ ( + )/2. O 3º passo é calcular calcular ( ) e compará-lo com ( ). ( ) ( )
24. São três as possibilidades: i. Se × = 0, a raiz já é . é, = . () ii iii ii. Senão, se × ( ) < 0, a raiz r está na metade da esquerda. Nesse caso o novo será (i. é, ↤ ). iii. Caso contrário, a raiz r está na metade da direita. Nesse caso ↤ ( passa a ser ).
25. No caso i, encerramos o processo, informando: A raiz é = . Nos outros dois casos, precisamos testar a precisão para decidir: A. Se | − | < e < , então colocamos ∗ ↤ e encerramos; B. Senão, retornamos ao 2º passo para mais uma iteração.
26. São duas condições de parada para anunciar ∗: • Proximidade da raiz: | − | < • Anulação aproximada do valor da função: < Um lembrete: não adianta exigir maior precisão que a “máquina” pode fornecer.
27. Escrevemos esta é a função para implementar o método da bisseção. Assinalei em vermelho as condições de parada.
28. O programa é para uma equação genérica () = 0. Os dados são: • Domínio de f • Expressão de f • Máximo de iterações • O limite para o erro Dados Chamada do método. Resultados Gráfico
29. Aqui testamos para achar 2, com a equação 2 − 2 = 0.
30. • Começamos com 0 = 1.4 e 1 = 1.5, sabendo que 2 está entre os dois • Portanto para = 1 − 0 , tanto | 0 − 2 | < como | 1 − 2 | < • Em seguida dividimos o intervalo à metade com 2 e garantimos que 2 está numa das metades • Portanto | 2 − 2 | < /2 • Continuando o processo, garantimos (prove!) que | +1 − 2 | < /2 Mas e a precisão na aproximação para 2 ? Então Loirinho (...), nosso procedimento foi o seguinte:
31. Além de dar uma estimativa para o erro, essa última desigualdade garante, também, que a sequência 0 , 1 , 2 , ⋯ , , ⋯ converge para 2, isto é lim →∞ = 2. Sim Mestra, pois de | +1 − 2 | < 1/2 garantimos que lim →∞ | +1 − 2 | < lim →∞ 1/2 = 0
32. Mestra, sabemos que 210 = 1024. Assim, | 10+1 − 2 | < 1 210 = 1 1024 < 0.001, Isto garante que 11 aproxima 2 com erro menor que um milésimo. Sim. E qual será o valor de k para o qual garantimos uma aproximação com a precisão do Single, de 6~7 casas decimais corretas?
33. Bem jovens, para uma precisão de 7 casas decimais, precisamos garantir que | +1 − 2 | < 0.5 × 10−7. Sim Mestra, e como | +1 − 2 | < 1/2, basta garantirmos que 1/2 < 0.5 × 10−7. Resolvendo essa inequação obtemos k = 25.
34. Como? Ora, de 1/2 < 0.5 × 10−7, segue 107 < 2−1 e portanto 7 log 10 < log(2−1). Então − 1 0.301030 > 7 e assim k = 25.
35. Agora vamos obter as raízes (se houver alguma) da equação 33 − 52 + 2 = 0 E como começamos Mestre? Não tenho a mínima ideia do lugar onde estão as raízes dessa equação!
36. Esse é um bom primeiro passo, meu Surfista! Bem Loirinha, uma ideia é subdividir um intervalo grande em pedacinhos e verificar se há variação de sinal nos pedacinhos. a b
37. Fiz uma função troca_sinal( ) para indicar os sub-intervalos:
38. Aqui estou chamando a troca_sinal( ) para a função = 33 − 52 + 2. Primeiro chamei com 50 sub-intervalos e depois repetí com 100. Minha conclusão é que ela só tem uma raiz real.
39. É, você só localizou uma raiz real. Entretanto, pela teoria, toda polinomial de grau 3 tem 3 raízes reais ou complexas. Pode ser que as outras duas sejam complexas conjugadas. Será que sua busca está certa?
40. Os resultados do método da bisseção, feito pelo Mestre, para a equação 33 − 52 + 2 = 0 e 4 casas após a vírgula, no intervalo [−1.6, −0.5].
41. Bem jovens, uma outra possibilidade é examinar o gráfico da função. O Professor fez um programa com o MatPlotLib que gera o gráfico de uma função num intervalo. Surfista, fiquei intrigada com a possibilidade de existirem mais raízes reais para aquela equação polinomial de grau 3. Como poderemos decidir?
42. Eis meu programa:
43. Observem o detalhe da utilização da função eval( ), uma função de Python, que transforma uma string num código executável.
44. Testando a função localizadora graf_func( ) com a função = 2 − 2, velha conhecida.
45. Agora vamos utilizar a função construída pelo Professor para localizar as raízes da equação 33 – 52 + 2 = 0 e tirar a dúvida da Loirinha!
46. Vou mostrar apenas uma vez como chamar a função graf_func( ) para desenhar o gráfico da função () = 33 – 52 + 2. Depois é igual.
47. Obtivemos este gráfico. Dele não concluo nada! Comece observando a diferença enorme entre as escalas horizontal e vertical.
48. Olhando para o gráfico, concluímos que não existem raízes fora do intervalo [-10. , 10.] É verdade Mestre. Vou traçar o gráfico nesse intervalo [-10. , 10.]
49. A diferença entre as escalas ainda é muito grande. Vejam! Reduza para o intervalo [-5. , 5.], pois não existem raízes fora dele, Loirinha
50. Reduzi, Mestre. E podemos reduzir ainda mais. Vemos que não há raízes fora de [-2. , 2.].
51. Ok Mestre. Agora podemos afirmar que as raízes estão no intervalo [-1. , 1.5]. Vou conferir!
52. Puxa! O Surfista localizou a raiz negativa, mas perdeu as duas positivas! Vou examinar o intervalo [0.5 , 1.5] .
53. Agora vemos que uma das raízes positivas provavelmente é 1.0. A outra está entre 1.2 e 1.3.
54. 3 ∗ (1 ∗∗ 3) – 5 ∗ (1 ∗∗ 2) + 2 = 3 – 5 + 2 = 0 Substituindo e fazendo as contas, conferimos que 1.0 é, de fato, raiz da equação 33 – 52 + 2 = 0. Veja só:
55. Para a 3ªraiz apliquei o método da bisseção feito pelo Mestre, entre 1.2 e 1.3 e obtive:
56. Pois é, só quando usei meu programa para 150 subintervalos é que localizei todas as trocas de sinal...
57. A Scipy oferece muitos métodos para achar raízes de equações. Vamos lá descobri-los.
58. Na API da SciPy temos:
59. Pule toda a parte referente à otimização indo direto para descoberta de raízes.
60. É o método da bisseção que eu já ensinei. A partir de agora, passem a usar esta implementação.
61. O restante dos parâmetros e o retorno.
62. Agora vamos estudar o método de Newton-Raphson.
63. Sempre que falam do Newton, lembro da derivada e da reta tangente. Sim Loirinha, ela “cola” na função perto do ponto de tangência! E o valor numérico da derivada é o coeficiente angular da reta tangente.
64. O método de Newton –Rhapson usa essa propriedade da reta tangente para gerar uma sequência 0 , 1 , 2 , … , , … de aproximações para a raiz.
65. ( ) f r raiz Veja como na figura, Loirinha. E, dado , como eu calculo +1 , Mestre? +1
66. ( ) f r Confira, Surfista, que equação da reta tangente por ( , ) é = ′( )( − ) + ( ) . E, como a Mestra mostrou, +1 = 0. +1
67. ( ) f r raiz Assim, fazendo = +1 na expressão da reta tangente, obtemos 0 = ′( )(+1 − ) + ( ). Então é só isolar +1 : +1 = − ( )/′( ), se ′( ) ≠ 0 +1
68. O processo iterativo é o seguinte: 1. Escolhemos o 1º valor 0 . É o “chute inicial”. Deve estar próximo da raiz; 2. Dado definimos +1 como o ponto onde a reta tangente por ( , ) corta o eixo-x; 3. Paramos quando a precisão for satisfatória; senão voltamos ao passo 2.
69. Portanto, o método de Newton-Raphson é definido por: ቐ +1 = − ′ , = 0,1, ⋯ 0 dado Na SciPy é tudo muito fácil, passamos: • A expressão da função ; • O chute inicial 0 ; • A expressão da derivada da função ; • A tolerância (erro) na aproximação; • Um número máximo de iterações (precaução).
70. O método de Newton na SciPy:
71. A continuação do Newton:
72. Um programinha para achar 2 pelo método de Newton. Vamos buscar raízes da equação 2 = 2. Assim = 2 − 2.
73. Então, Surfista, faça um programa para achar raízes de uma equação = 0 pelo método de Newton, sem usar a SciPy, que imprima cada uma das iterações. É facílimo Mestra !
74. Vou retornar à definição de convergência de sequências. Como já disse, precisaremos trabalhar com “épsilons e deltas”.
75. Na determinação de raízes de uma equação através de métodos iterativos temos seis entidades envolvidas: 1. O método utilizado, 2. A equação, dada por () = 0, 3. A função : , → ℝ, 4. A que estamos procurando, 5. A sequência gerada pelo método que estamos usando, tal que → , 6. A sequência auxiliar definida pelos valores = ( ).
76. Por exemplo: • Bisseção, • Newton-Raphson, • Secante Métodos 2 − 2 = 0 Equação ↦ 2 − 2 Função 2, − 2 Raízes | Sequência | = ( ) Sequência
77. Em todos os métodos, a função f é contínua e funções contínuas “empurram a convergência” no sentido que → ⟹ → (). Todos os métodos iterativos para determinação de uma de uma equação () = 0 envolvem a geração de uma sequência ( ) tal que → .
78. Como = 0 temos: → ⇒ = ( ) → 0 | = ( ) | → = ( ) → 0
79. Por garantia, conferimos também se = − 0 < | = ( ) Portanto, além de testarmos se − < |
80. Isto ficou claro. O que não entendi, é a outra condição de parada. Não temos que testar se | − | < , conforme fizemos no caso de = 2, | − 2 | < ? Pois é, minha filha, esse exemplo foi didático. Já sabíamos que = 2 é a raiz. No caso geral não – a raiz é o que buscamos!
81. E como procedemos? Como não sabemos quem é r, não temos como testar a proximidade, | − | < ,
82. Pois é Loirinha, este é um problema crucial, cuja resposta foi dada por Cauchy.
83. Uma sequência ( ) é de Cauchy (uma homenagem a ele) quando, e apenas quando, satisfaz a propriedade: Não estou vendo conexão entre esse novo simbolês e nosso problema! ∀ > 0, ∃ ∈ ℕ tal que , ≥ ⟹ | – | < .
84. r ( ) 0 1 2 Bem Surfista, é só trocar esta fala, já dita pelo Filósofo, pela do Sherlock, na próxima transparência. Repetindo: A sequência todinha, exceto alguns termos, está aqui!
85. ( ) 0 1 2 A sequência todinha, exceto uma quantidade finita de termos, entra numa caixinha de largura ε . E não sai mais! Claro! Então a raiz r também estará nessa caixinha!
86. Augustin-Louis Cauchy 21/08/1789 23/05/1857 Bem, eu aprendi tudo isto lendo o “Cours d’analyse”, de Cauchy, publicado em 1821...
87. Pois é Surfista, mas a sua conclusão só vale naqueles mundos que os matemáticos chamam de completos. Mas tranquilize-se, meu jovem, os espaços euclidianos (classe da qual o mundo em que vivemos faz parte) são completos.
88. Neles, podemos conferir a convergência comparando a proximidade de dois termos quaisquer, isto é, testando se | – | < , para cada valor escolhido de . Nos espaços completos, em particular nos euclideanos, toda sequência de Cauchy é convergente e vice-versa.
89. Em particular, testar se +1 – < é uma forma de garantir a proximidade de +1 com a raiz. Juro que não entendi! Ah Colega, é só fazer = + 1 e = .
90. Aristóteles 384-322 410 390 370 350 310 290 1915 430 330 1910 1890 1930 1950 Einstein 1879-1955 Platão 428-347 Euclides 360-295 ±300 Muito bem lembrado, Filósofo. Mas até ele publicar sua da Teoria da Relatividade Geral (1915 dC), todos acreditavam nisso. Mestra, não posso deixar de lembrar que Einstein provou que o espaço em que vivemos não é euclidiano.
91. Não existe uma prova única de convergência, válida para todos os métodos iterativos. É necessário efetuar uma prova específica para cada método. Em particular, a prova de convergência do método de Newton-Rhapson pode ser encontrada em diversos livros de Cálculo/Análise Numérica.
92. Logicamente NÃO, Loirinha. A convergência só é garantida quando você consegue exibir um valor de N tal que , > ⇒ – < , para cada valor escolhido da precisão . Não apenas para = 0.5 × 10−7. Mas Mestre, quando obtenho numericamente, que | +1 − | < 0.5 × 10−7, para algum valor de k, já não garanti a convergência?
93. Para cada caixinha que a Mestra escolher, por minúscula que seja, você sempre conseguirá guardar dentro dela a sequência todinha, exceto alguns termos. Repita a frase do Mestre de forma mais saborosa, Sherlock.
94. 1. Ao perceber um comportamento indicativo de convergência, acreditamos que a sequência é convergente (logo é Cauchy). 2. Para cada valor de n calculamos a diferença |+1 − |. 3. Paramos quando conseguimos +1 − < e | +1 | < para algum valor escolhido da precisão . 4. Assumimos +1 como aproximação ∗ para r. Na realidade, na prática, fazemos o seguinte: O passo 1 da lista não é matemática, apenas intuição. O passo 3 é passível de contra-exemplos ...
95. Essa receita do Mestre não é respeitável. Quem diria! Mas é o que há, meu jovem. A própria Ciência não se garante com o Princípio da Indução. Procure descobrir o que falaram sobre o tema filósofos como David Hume, Thomas Kuhn, Karl Popper e outros.
96. Tchau, até a próxima aula!