Esta é uma pré-visualização de arquivo. Entre para ver o arquivo original
01_Introducao_v2_Folheto.pdf 06/08/2019 1 Cálculo Numérico Aula Inicial Prof. Marconi de Arruda Pereira marconi@ufsj.edu.br Dados Gerais Professor: Marconi de Arruda Pereira marconi@ufsj.edu.br Sala 101 - Bloco 3 Horário de atendimento: Terças das 13:15 às 16:15 Dados Gerais Objetivo Introduzir o aluno na área do Cálculo Numérico e da Análise Numérica, tornando-o capaz de analisar e aplicar algoritmos numéricos em problemas reais, codificando-os em uma linguagem de alto nível a fim de resolver problemas de pequeno e médio porte em Ciência e Tecnologia. Dados Gerais Ementa Teoria de erros. Zeros de funções e zeros reais de polinômios. Solução de sistemas lineares: métodos diretos e iterativos. Ajuste de curvas. Interpolação. Integração numérica. Resolução numérica de equações diferenciais ordinárias. Exemplos de aplicações do Cálculo Numérico na Engenharia. Aulas práticas em laboratório. Dados Gerais Conteúdo 1. Introdução 2. Teoria de erros Representação de números Erros absolutos e relativos Erros de arredondamento e truncamento Análise de erros Dados Gerais Conteúdo 3. Zeros de Funções Isolamento das raízes Critério de parada Métodos iterativos de cálculo de raízes (Bisseção, Posição Falsa, Ponto Fixo, Newton e Secante) Comparação entre os métodos 4. Solução de sistemas lineares Métodos Diretos (Eliminação de Gauss e Decomposição LU) Métodos Iterativos (Gauss-Jacobi e Gauss-Siedel) 06/08/2019 2 Dados Gerais Conteúdo 5. Interpolação Métodos de Interpolação (Linear , Quadrática , Lagrange, Newton, Gregory-Newton) Estudo do Erro 6. Ajuste de Curvas Método dos Quadrados Mínimos (Caso Linear e Não- Linear) Dados Gerais Conteúdo 7. Integração numérica Integração simples (Regra dos Trapézios, Primeira Regra de Simpson, Segunda Regra de Simpson) Integração dupla 8. Resolução numérica de equações diferenciais ordinárias Problemas de Valor Inicial Método de Euler Dados Gerais Bibliografia Básica 1. RUGGIERO, Márcia A. G., LOPES, Vera L. R. Cálculo Numérico - Aspectos teóricos e computacionais. 2a Ed. Pearson, 1996; 2. FRANCO, Neide Bertoldi. Cálculo Numérico. 1a Ed. Prentice Hall, 2006; 3. CHAPRA, Steven C., CANALE, Raymond P. Métodos Numéricos para a Engenharia. 5ª Ed. MCGRAW-HILL BRASIL, 2008. Dados Gerais Bibliografia Complementar 1. CAMPOS, filho, Frederico F. Algoritmos Numéricos, 2.ed., Rio de Janeiro: LTC, 2007; 2. BARROSO, Leônidas, BARROSO, Magali Maria de Araújo, CAMPOS FILHO, Frederico Ferreira. Cálculo Numérico com Aplicações. 2a Ed. Harbra, 1987; 3. SPERANDIO, Décio, MENDES, João T., SILVA, Luiz H. M. Cálculo numérico - características matemáticas e computacionais dos métodos numéricos. 1a Ed. Prentice Hall. 2003. O que é combinado... Chamada Serão realizadas em todas as aulas A presença e principalmente a participação do aluno é fundamental Exercícios práticos Acontecerão, idealmente, todas as semanas nas quais houver conteúdo expositivo Listas são verificadas Penalidade de 20% por atraso de até 7 dias (máximo 2 semanas de atraso) Atividade Extra: Clube de leitura http://petdpcfc.wixsite.com/ufsj/grupo-de-leitura O que é combinado... Exercícios Não copie, Não compre Compra Nota zero Cópia Origem e Cópia recebem a mesma nota, zero As normas acadêmicas prevêem reprovação e advertência Provas Uso de calculadora científica é permitido. HP, não. Estude agora para não precisar depois. 06/08/2019 3 Avaliações Exercícios Semanais (2,5 total) Resultado final em 04/12 4 Provas (2,5 pontos) 18/09 23/10 20/11 04/12 -> Substitutiva (substitui prova ou exercícios) Segunda chamada (em acordo com os Procedimentos Acadêmicos) 11/12 Matlab/Octave Poderosas ferramentas matemáticas que nos auxiliarão no aprendizado do Cálculo Numérico Laboratórios Trabalhos Práticos Disponibiliza funções numéricas Possibilita a programação de novas funções Cálculo Numérico Introdução Prof. Marconi de Arruda Pereira marconi@ufsj.edu.br material produzido pelo prof. Wellington Passos O que é Cálculo Numérico? O Cálculo Numérico corresponde a um conjunto de ferramentas ou métodos usados para se obter a solução de problemas matemáticos de forma aproximada, mas com um grau crescente de exatidão Esses métodos se aplicam principalmente a problemas que não apresentam uma solução exata, portanto precisam ser resolvidos numericamente Por que utilizar? Um problema de Matemática pode ser resolvido analiticamente, mas esse método pode se tornar impraticável com o aumento do tamanho do problema Exemplo: solução de sistemas de equações lineares (soluções químicas, redes elétricas etc.) Por que utilizar? O problema não tem solução analítica Exemplos: a) não tem primitiva em forma simples; b) Equações diferenciais parciais não lineares podem ser resolvidas analiticamente só em casos particulares. 06/08/2019 4 Solucionar problemas técnicos através de métodos numéricos, usando um modelo matemático Função do Calculo Numérico na Engenharia Exemplo de aplicação Calcular tensões dos nós do circuito elétrico No nó 1, pela lei de Kirchhoff Exemplo de aplicação O problema é resolvido a partir de um sistema linear de quatro equações e quatro variáveis V1, V2, V3 e V4. 0 254 0 0 4 3 2 1 3201 61323 0143 1126 V V V V Resolução de problemas No exemplo anterior Problema real: determinar tensões nos nós dos circuitos Levantamento de dados: valores das resistências e tensões nos pontos A e B Construir modelo matemático: montar equações e criar as matrizes a partir delas Escolher método numérico: Decomposição LU, Decomposição de Cholesky, Fatoração LDLT, Método de Jacobi, etc Implementar Método: criar e processar programa Analisar resultados e verificar se o modelo matemático ou o método numérico precisam ser alterados. Influência dos erros nas soluções Exemplo 1: Falha no lançamento de mísseis (25/02/1991 – Guerra do Golfo – míssil Patriot) Limitação na representação numérica (24 bits) Erro de 0,34 s no cálculo do tempo de lançamento Resultado: 28 mortos e aprox. 100 feridos 06/08/2019 5 Influência dos Erros nas Soluções Exemplo 2: Explosão de foguetes (04/06/1996 – Guiana Francesa – foguete Ariane 5) Limitação na representação numérica (64 bits/ 16 bits) Erro de trajetória levou a explosão 36,7 s após o lançamento Prejuízo = US$ 500 milhões em equipamentos Outros exemplos Desastres Causados por erros de cálculo numérico http://www.ima.umn.edu/~arnold/disasters/ 02_Erros_Folheto.pdf 18/02/2019 1 Cálculo Numérico Teoria dos Erros Autor: Prof. Wellington Passos de Paula Adptado por: Marconi de Arruda Pereira Programa 1. Conceitos Básicos a) Representação de números b) Conversão de números c) Aritmética de ponto flutuante 2. Erros a) Erros absolutos e relativos b) Erros de arredondamento e truncamento c) Análise de erros Cálculo Numérico Teoria dos Erros – Conceitos Básicos Autor: Prof. Wellington Passos de Paula Adptado por: Marconi de Arruda Pereira Representação de números Sistema Decimal (10) 10 dígitos disponíveis [0,1,2, ... ,9] “Posição” indica potência positiva de 10 5432 = 5x103 + 4x102 + 3x101 + 2x100 Sistema Binário (2) 2 “bits” disponíveis [0,1] “Posição” indica potência positiva de 2 1011 na base 2 = 1x23 + 0x22 + 1x21 + 1x20 8+0+2+1 = 11 na base decimal Representação de números Fórmula Geral Base : Logo, a decomposição polinomial do número é dada por: Exemplo: Dado , temos que: , k j,..., 0 ,)...( 0121 aaaaa jj - )1(0 - ka 0 0 1 1 2 2 1 1 ... ´´´´´ - - aaaaa j j j j ,)...( 0121 aaaaa jj - 0123 1091041081066849 ´´´´ 10 Representação de números Representação Números Fracionários Base Decimal (10) “Posição” da parte inteira indica potência positiva de 10 Potência negativa de 10 para parte fracionária 54,32 = 5x101 + 4x100 + 3x10-1 + 2x10-2 Base Binária (2) “Posição” da parte inteira indica potência positiva de 2 Potência negativa de 2 para parte fracionária 10,11 na base 2 = 1x21 + 0x20 + 1x2-1 + 1x2-2 2+0+1/2+1/4 = 2,75 na base decimal 18/02/2019 2 Outros sistemas de numeração Maior interesse em decimal (10) Nossa anatomia e cultura Binário (2) – uso nos computadores Outros Sistemas Octal (8), {0,1,2, ... , 7} Hexadecimal (16), {0,1,2, ... , 9, A,B,C,D,E,F} Dodecimal (relógio, calendário) Alguns sistemas numéricos Conversão de números - inteiros Binário para decimal Já visto (1011)2 = 1x2 3 + 0x22 + 1x21 + 1x20 = (11)10 Inteiro decimal para binário Divisão inteira (do quociente) sucessiva por 2, até que este seja = 0 ou 1 Binário = composição do último quociente (Bit Mais Significativo – BMS) com os restos das divisões (primeiro resto é bit menos significativo – bms) Em inglês, Most Significant Bit – MSB e least significat bit – lsb, respectivamente. Conversão de números - inteiros Exemplo: Converter 30 decimal para binário Binário = BMS ... bms = 1 1 1 1 0 1 1 1 1 0 = 1x24 + 1x23 + 1x22 + 1x21 + 0x20 = 16 + 8 + 4 + 2 + 0 = 30 decimal Conversão de inteiros entre sistemas Procedimentos Básicos Decimal Binário - Divisões sucessivas pela base do sistema para o qual se deseja converter o número Binário Decimal - Decomposição polinomial do número a ser convertido Conversão de inteiros entre sistemas 18/02/2019 3 Conversão de fração Base 2 para Base 10 Já visto (Decomposição Polinomial) (10,101)2 = 1x2 1 + 0x20 + 1x2-1 + 0x2-2 + 1x2-3 = = 2 + 0 + 1/2 + 0 + 1/8 = (2,625)10 Conversão de fração Base 10 para Base 2 Deve-se multiplicar parte fracionária por 2 até que parte fracionária do resultado seja 0 (zero) X,XXX Bits da parte fracionária do número binário são obtidos das partes inteiras geradas após as multiplicações do número fracionário na base 10 X,XXX Bit imediatamente à direita da vírgula = Parte inteira da primeira multiplicação Não há inversão na ordem dos bits encontrados Conversão de fração Exemplo: converter 0,625 decimal para binário 0,625 x 2 = 1,25, logo a primeira casa fracionária do número binário é 1; nova fração (resto) é 0,25 (agregamos o bit 1 ao número na base 2) 0,25 x 2 = 0,5 segunda casa do número binário é 0; nova fração (resto) é 0,5 (pois já agregamos o bit 0 ao numero na base 2) 0,5 x 2 = 1,0 terceira casa é 1; nova fração (resto) é 0,0 (pois já agregamos o bit 1 ao numero na base 2) Resultado: 0,62510 = 0,1012 Conversão partes inteira e fracionária juntas Para converter um número com parte inteira e parte fracionária, faça a conversão de cada parte separadamente Conversão partes inteira e fracionária juntas (8,375)10 = ( ? )2 Exercícios Transforme em binário: 5,8 Resposta: 5,8 = 101,11001100... , uma dízima. 11,6 Resposta: 11,6 = 1011,10011001100... a vírgula foi deslocada uma casa para a direita, pois 11,6 = 2 x 5,8 18/02/2019 4 Aritmética de ponto flutuante Representação pode variar (“flutuar”) a posição da vírgula, ajustando potência da base. 54,32 = 54,32 x 100 = 5,432 x 101 = 0,5432 x 102 = 5432,0 x 10-2 Forma normalizada utiliza um único dígito antes da vírgula ( 0 ), e garante que o primeiro dígito depois da vírgula seja diferente de 0 Exemplo: 0,5432 x 102 Aritmética de ponto flutuante No sistema binário: 11010 = 11,010 x 23 = 0,11010 x 25 = 0,0011010 x 27 No caso dos números serem armazenados em um computador, os expoentes serão também gravados na base dois 11,010 x (2)11 = 0,11010 x (2)101 = 0,0011010 x (2)111 Na representação normalizada, há apenas um dígito antes da vírgula ( 0 ) Exemplo: 0,11010 x (2)101 Aritmética de ponto flutuante Algumas definições No número 0,11010 x (2)101 , tomado como referência: 0,11010 = significando (ou “mantissa”) 101 = expoente Observações A base binária não precisa ser explicitada (o computador usa sempre a mesma) O “0” antes da vírgula, na representação normalizada – se esta for adotada, também pode ficar implícito, economizando um bit (“bit escondido”). Representação aritmética de ponto flutuante no computador onde: é a base em que o computador opera; é o número de dígitos na mantissa é o expoente (inteiro com sinal) e tddd ´ )...(. 21 t ;01 d,,...,1),1(0 tjd j - e Representação aritmética de ponto flutuante no computador O número de bits disponíveis para representar os números no computador não é infinito O padrão IEEE 754 para ponto (vírgula) flutuante é a representação mais comum para números reais em computadores de hoje, incluindo PC's compatíveis com Intel, Macintosh, e a maioria das plataformas Unix/Linux. O padrão (ou norma) IEEE 754 define dois formatos básicos para os números em ponto flutuante: O formato ou precisão simples, com 32 bits; e, O duplo com 64 bits Sinal: 0 = + e 1 = - Combinações: Sinal + Expoente + Mantissa Padrão IEEE 754 para floats Sinal Expoente(+/-) Mantissa Simples (32bits) 1 [bit31] 8 [bits30-23] 23 [bits22-00] Dupla (64 bits) 1 [bit63] 11 [bits62-52] 52 [bits51-00] 18/02/2019 5 Limitações na representação de floats A quantidade finita de bits na representação pode implicar nos seguintes erros: Truncamento ou Arredondamento Overflow Underflow Limitações na representação de floats Exemplo: Máquina no seguinte sistema: Logo o formato dos números nesse sistema: Menor valor representado em módulo: Maior valor representado em módulo: .5,5;3;10 - et 5,5,0,90,10.0 1321 -´ eddddd j e 65 1010100.0 -- ´m 9990010999.0 5 ´M Limitações na representação de floats Situações possíveis: a) . Número contém 5 dígitos na mantissa Possíveis Soluções: Truncamento: Arredondamento: Assunto do próximo tópico 31023589.089.235 ´x 310235.0 ´ 310236.0 ´ Limitações na representação de floats Situações possíveis: b) . Expoente não pode ser representado na máquina pois é menor que o mínimo (-5) Erro de underflow c) . Expoente não pode ser representado na máquina pois é maior que o máximo (5) Erro de overflow 710345.0 -´x 910875.0 ´x Limitações na representação de floats Considere ]4,4[;3;10 - et x arredondamento truncamento 1.25 10.053 -253.15 2.71828 0.000002 Underflow Expoente<-4 817235.89 Overflow Expoente>+4 110125.0 ´ 110125.0 ´ 210100.0 ´210101.0 ´ 310253.0 ´- 310253.0 ´- 110272.0 ´ 110271.0 ´ Exercícios Considere uma máquina com sistema de representação de números definido por: base 10, precisão de 4 dígitos na mantissa e expoente no intervalo: [-6; 6]. Pede-se: a) Qual o menor e o maior número em módulo representado nesta máquina? Menor: 0.1000x10-6 = 10-7, Maior: 0.9999x106 = 999900 b) Como será representado o número 189,27 nesta máquina se for usado o arredondamento? E se for usado o truncamento? Trunc.: 0.1892x103, Arred.: 0.1893x103 c) Se a = 2578 e b = 0,6 qual o resultado de a + b se for usado o arredondamento? E se for usado o truncamento? Trunc.: 0.2578x104, Arred.: 0.2579x104 18/02/2019 6 Cálculo Numérico Teoria dos Erros – Erros Autor: Prof. Wellington Passos de Paula Adaptado por: Marconi de Arruda Pereira Erros - Tipos Precisão Absoluto Relativo Representação Arredondamento Truncamento Erro Absoluto Diferença entre o valor exato de um número e o seu valor aproximado (em módulo) |x|xEAx - EAx só poderá ser determinado se x for conhecido com exatidão Na prática, costuma-se trabalhar com um limitante superior para o erro, ao invés do próprio erro ( |E | < ε, sendo ε é o limitante) Ex: Para (3,14; 3,15) = 3,1415926535897932384626433832795028841971693… EAπ = π- π < 0,01 Erro Absoluto - Considerações Erro Absoluto - Considerações Ex.: Sejam a = 3876,373 e e = 1,373 Considerando-se a parte inteira de a como ā o erro absoluto será: e a parte inteira de e, ē, o erro absoluto será: 0,373aaEAa - 0,373eeEAe - Erro Absoluto - Considerações Obviamente, o resultado do erro absoluto é o mesmo nos dois casos Podemos então dizer que a e e estão representados com a mesma precisão? Não, pois o peso da aproximação em e é maior do que em a Erro absoluto não é suficiente para descrever a precisão de um cálculo 18/02/2019 7 Erro Relativo Razão entre o erro absoluto e o valor aproximado do número considerado (em módulo) |x| EA |x| |x|x ER xx - O erro relativo pode ser considerada como uma estimativa da precisão do número: ERx x 100 = Erro Percentual O valor exato da precisão = |Erro absoluto|/|valor exato| 59 -´ 100,000096 3876 0,373 ERa 140 -´ 10,373 1 0,373 ERe Erro Relativo - Considerações Ex. : Cálculo do erro relativo na representação dos números ā = 2112,9 e ē = 5,3, sendo |EA| < 0,1 Conclusão: a é representado com maior precisão do que e Erro Relativo - Considerações 5107,4 9,2112 1,0 -´ - a aa ERa 02,0 3,5 1,0 - e ee ERe Erros de Arredondamento Ex. Cálculo de utilizando uma calculadora digital: Valor apresentado: 1,4142136 Valor real: 1,41421356... Não há forma de representação digital (por máquinas) de números irracionais com uma quantidade finita de algarismos. A calculadora apresenta uma aproximação do número. Erro de Arredondamento. 2 Erros de Truncamento É o descarte dos dígitos finais de uma representação exata por limitações de representação em vírgula flutuante: Ex.: Representação truncada de em vírgula flutuante com 7 dígitos Valor apresentado: 1,4142135 Valor real: 1,41421356... 2 Representação aritmética de ponto flutuante no computador – Relembrando... onde: é a base em que o computador opera; é o número de dígitos na mantissa é o expoente (inteiro com sinal) e tddd ´ )...(. 21 t ;01 d,,...,1),1(0 tjd j - e 18/02/2019 8 Para t = 4, e = 3 e x = 234,57: x = 0,2345 x 103 + 0,7 x 10-1 fx = 0,2345 gx = 0,7 Erros de Truncamento e Arredondamento em um sistema de aritmética de ponto flutuante: Em um sistema que opera em ponto flutuante de t dígitos na base 10, e seja x: x = fx x 10 e + gx x 10 e-t (0,1 fx 1 e 0,1 gx 1) Erros de Arredondamento e Truncamento Para t = 5, e = 4 e x = 1234,568 x = 0,12345 x 104 + 0,68 x 10-1 fx = 0,12345 gx = 0,68 Erros - Truncamento No truncamento, gx x 10 e-t é desprezado e e , visto que |gx| < 1 pois 0,1 é o menor valor possível para fx 1t e te e te e x te xx x 10 1010 10 100,1 10 10f 10g x EA ER - - --- ´´ ´ ´ ´ 10,1 e x 10fx ´ tete xx 1010gEA -- ´ te x e x 10g10fx -´´ e x te x e xx f10gfxxEA 1010 ´-´´- - No arredondamento simétrico (forma mais utilizada): , se (gx é desprezado) , se (soma 1 ao último dígito de fx) Erros – Arredondamento ´ ´ -tee x e x 1010f 10f x 2 1 g x 2 1 g x Erros - Arredondamento Se , então: , visto que |gx| < 1/2 1t e te e te e x te xx x 10 2 1 10 10 100,1 100,5 10f 10g x EA ER - - --- ´ ´´ ´ ´ ´ ´ ´ 1100,1 2/1 tete xx 10 2 1 10gEA -- ´´ 2 1 g x e x te x e xx f10gfxxEA 1010 ´-´´- - Erros – Arredondamento Se , então: e e x te tee x te x x 10f 101/2 1010f 101/2 x EA ER ´ ´ ´ ´ - - - tetex tete xx 10 2 1 101g1010gEA ---- ´´--´ 2 1 g x teextexexx 1010f10g10fxxEA -- ´-´´- Erros de Truncamento e Arredondamento em um sistema de aritmética de ponto flutuante: Sistema operando em ponto flutuante - Base 10 Erro de Truncamento e Erro de Arredondamento e Arredondamento gera erros menores, mas aumenta o tempo de execução uso do Truncamento te x 10EA - 1t x 10ER - Arredondamento e Truncamento 18/02/2019 9 Sistema de aritmética de ponto flutuante de 4 dígitos, precisão dupla Ex.: Seja x = 0,937 x104 e y = 0,1272 x102. Calcular x+y. Alinhamento dos pontos decimais antes da soma ( Alinhar sempre para o maior expoente dentre os operadores ) x = 0,937 x 104 e y = 0,001272 x 104, x+y = 0,937 x 104 + 0,001272 x 104, x+y = 0,938272 x 104 Resultado com 4 dígitos Arredondamento: x+y = 0,9383 x 104 Truncamento: x+y = 0,9382 x 104 Análise de Erros Sistema de aritmética de ponto flutuante de 4 dígitos, precisão dupla Ex. : Seja x = 0,937 x 104 e y = 0,1272 x102. Calcular a multiplicação x*y x*y = (0,937 x 104) x (0,1272 x 102) x*y = (0,937 x 0,1272) x 106 x*y = 0,1191864 x 106 Resultado com 4 dígitos Arredondamento: x*y = 0,1192 x 106 Truncamento: x*y = 0,1191 x 106 Análise de Erros Análise de Erros Considerações Ainda que as parcelas ou fatores de uma operação possam ser representados exatamente no sistema, não se pode esperar que o resultado armazenado seja exato. x e y tinham representação exata, mas os resultados x+y e x.y tiveram representação aproximada. Durante as operações aritméticas de um método, os erros dos operandos produzem um erro no resultado da operação Propagação ao longo do processo Análise de Erros – Propagação Ex. : Sejam as operações a seguir processadas em uma máquina com 4 dígitos significativos e, sendo a = 0,3491 x 104 e b = 0,2345 x 100, qual o valor da expressão: (b + a) − a = b + (a − a) ? (b + a) − a = (0,2345 x100+0,3491x104) − 0,3491x104 = (0,00002345 x104+0,3491x104) − 0,3491x104 (0,34912345 x104) − 0,3491x104 (arredodamento) 0,3491 x 104 − 0,3491 x104 = 0,0000 b + (a − a) = 0,2345x100 + (0,3491 x 104 −0,3491x104)= 0,2345 x 100 +(0,0000 x 104)= 0,2345 x 100 Análise de Erros – Propagação Os dois resultados são diferentes, quando não deveriam ser. (b + a) − a = 0,0000 e b + (a − a) = 0,2345 x 100 Causa Arredondamento da adição (b + a), a qual tem 8 dígitos A máquina só armazena 4 dígitos (desprezando os menos significativos) Análise de Erros – Propagação É preciso atenção na resolução numérica de um problema. É Importante o conhecimento dos efeitos da propagação de erros. Determinação do erro final de uma operação; Conhecimento da sensibilidade de um determinado problema ou método numérico. 18/02/2019 10 Análise de Erros – Propagação Análise dos Erros Absoluto e Relativo. Existem expressões para o determinação dos erros nas operações aritméticas. Erros podem ocorrer na representação das parcelas ou fatores, assim como no resultado da operação. Supondo um erro final arredondado, sendo x e y, tais que: yx EAy yEAxx e Análise de Erros – Propagação Adição Erro Absoluto Erro Relativo yx y ER yx x ER yx y y EA yx x x EA ER yx yx yx )EA(EA)yx( )EAy() EAx(yx yx yx yxyx EAEAEA yx EA yx EA yx EAEA yx EA ER yxyxyxyx Análise de Erros – Propagação Subtração Erro Absoluto Erro Relativo )EA(EA)yx( )EAy()EAx(yx yx yx -- -- yxyx EAEAEA -- - - - - - - - yx y ER yx x ER yx y y EA yx x x EA ER yx yx yx - - - - - - - - yx EA yx EA yx EAEA yx EA ER yxyxyx yx Análise de Erros – Propagação Multiplicação Erro Absoluto Erro Relativo muito pequeno yxxyyx EAEAEAyEAxyxEAyEAxx.y ´´´´´ xyyx EAyEAxyxEAyEAxx.y ´´´´ yxx.y ERERER xyx.y EAyEAxEA y EA x EA xy EAy xy EAx xy EAyEAx xy EA ER yxxyxyyx x.y . Análise de Erros – Propagação Divisão Erro Absoluto ´ y EA 1 1 y EAx EAy EAx y x y x y x yyy EAy y y y EAyy y EA 1 1 y ´ ´ ´ 1111 Simplificação:: (desprezam-se os termos de potência >1) ´ y EA 1 1 y EAx EAy EAx y x y x y x ... y EA y EA y EA 1 y EA 1 1 3 y 2 yy y - - -´ -´ y EA y EA y x y EA y EAx y x yxyx 11 Análise de Erros – Propagação Divisão Erro Absoluto -´ -´ y EA y EA y x y EA y EAx y x yxyx 11 2 yx y EAx y EA y x y x - 2y EAEA y EA y EAx y x y x yxx 2 y ´ -- muito pequeno 2 yx yx y EAxEAy EA - / 18/02/2019 11 Análise de Erros – Propagação Divisão Erro Relativo - ´ x y y EAxEAy x y EA y x EA ER 2 yx yx yx x/y / / yx yx x/y ERER y EA x EA ER -- Análise de Erros - Propagação Erro Relativo da Adição É a soma dos erros relativos de cada parcela, ponderados pela participação de cada parcela no total da soma. Erro Relativo da Subtração É a diferença entre os erros relativos do minuendo e do subtraendo, ponderados pela participação de cada parcela no resultado da subtração. yx y ER yx x ERER yxyx - - - - yx y ER yx x ERER yxyx Análise de Erros - Propagação Erro Relativo da Multiplicação Soma dos erros relativos dos fatores. Erro Relativo da Divisão Diferença entre os erros relativos do dividendo e do divisor yxx.y ERERER yxx/y ERERER - Análise de Erros - Propagação Nos erros anteriormente formulados, ainda consideramos o erro de arredondamento ou truncamento no resultado final A análise completa da propagação do erro se faz considerando os erros nas parcelas ou fatores e no resultado de cada operação efetuada Ex.: Dada a soma x+y (x e y representados exatamente), faça o cálculo de ER(x+y) Como x e y são exatamente representados, ERx+y se resume ao Erro Relativo de Arredondamento (RA) no resultado da soma. EAx= EAy = 0, EAx+y = 0 1t yx 10 2 1 RAER - ´ RA yx EA ER yx yx Análise de Erros - Propagação RAER yx Análise de Erros - Propagação Sistema de aritmética de ponto flutuante de 4 dígitos, precisão dupla Ex.: Seja x = 0,937 x104, y = 0,1272 x102 e z = 0,231 x101, calcular x+y+z e ER(x+y+z), sabendo que x, y e z estão exatamente representados. Solução: Alinhando as vírgulas decimais ( Alinhar sempre para o maior expoente dentre os operadores ) : x = 0,937000 x104 y = 0,001272 x104 e z = 0,000231 x104 18/02/2019 12 Análise de Erros - Propagação Ex.: Seja x = 0,937 x104, y = 0,1272 x102 e z = 0,231 x 101, calcular x+y+z e ER(x+y+z), sabendo que x, y e z estão exatamente representados. Solução: A soma é feita por partes: (x+y)+z x+y = 0,937000 x104 + 0,001272 x104 x+y = 0,938272 x104 (arredondamento) x+y = 0,9383 x 104 = s s+z = 0,9383 x 104 + 0,000231 x 104 s+z = 0,938531 x 104 (arredondamento) x+y+z = 0,9385 x 104 Análise de Erros - Propagação Solução: s = x+y = então s = x + y = 0,9383 x 104 Cálculo do Erro Relativo: EAx=EAy=0, ERx+y=0syxs RAyx y ER yx x ERER ss RAER Análise de Erros - Propagação Solução: EAz=0, ERz=0 RA zyx z ER zyx yx ERER zszyx RA zyx yx ERER szyx 1 zyx yx RA RA zs z ER zs s ERER zszyx RA zyx yx RAER szyx Análise de Erros - Propagação Solução: 3 zyx 0,9998.10ER - 1t zyx 10 2 1 1 zyx yx ER - ´ 3 4 4 109385,0 109383,0 - ´ ´ ´ 10 2 1 1ER zyx 1 zyx yx RARA zyx yx RAER szyx Análise de Erros - Propagação Ex. : Supondo que u é representado em um computador por ū, que é obtido por arredondamento. Obter os limites superiores para os erros relativos de v = 2ū e w = ū + ū. Análise de Erros - Propagação Ex. : Solução: 1t u 10ER - 2 uv 2 RAERERER uu 22 1 2 10 2 1 2 -´´ t u ER RARARA 2 18/02/2019 13 Ex. : Solução: Análise de Erros - Propagação 1t vw 10ERER - uuw RA uu u ER uu u ERER uuw 11 1010 2 1 22 -- ´´´ ttw RAER RA uu u RAERw ´ 2 ´ RA u u RAERw 2 2 RA´2 Exercício Considere uma máquina cujo sistema de representação de números é definido por . Tal máquina utiliza o arredondamento para os dígitos na mantissa. Os números x = 8543 e y = 2477 foram utilizados em algumas operações nesta máquina. Assim, faça o que se pede: a) Calcule os erros absolutos (EA) e erros relativos (ER) envolvidos no processo de utilização da máquina para cada número x e y. Resposta: et 3,10 ]5,5[-e x 0,854´104 EAx 0, 0003´10 4 ERx 3, 513´10 -4 y 0,248´104 EAy 0, 0007´10 4 ERy 1, 210´10 -3 Exercício Considere uma máquina cujo sistema de representação de números é definido por . Tal máquina utiliza o arredondamento para os dígitos na mantissa. Os números x = 8543 e y = 2477 foram utilizados em algumas operações nesta máquina. Assim, faça o que se pede: b) Após a realização das operações x+y e x*y, foi percebido que uma das duas operações resultava no erro relativo maior. Qual foi? Resposta: Erro da multiplicação é maior et 3,10 ]5,5[-e RAER yx ´ - 410445,5 RAER yx ´ - ´ 410613,15 02_Erros_OBS-NumerosDiscretos.pdf 210 31 2 3 4 3 8 1 4 3 2 03_Raizes_Folheto.pdf 18/02/2019 1 Zeros Reais de Funções Reais Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Programa 1. Introdução 2. Isolamento das raízes 3. Refinamento a) Critério de parada b) Métodos iterativos c) Comparação entre os métodos Zeros Reais de Funções Reais – Introdução Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Zeros de funções reais - Objetivos Estudar métodos numéricos para a resolução de equações não lineares (determinar a(s) raiz(es) de uma função f(x), ou seja, encontrar o(s) valor(es) de x tal que f(x) = 0) Fundamentar a necessidade de uso de métodos numéricos para a resolução de equações não lineares Discutir o princípio básico que rege os métodos numéricos para a resolução de equações não lineares Apresentar e discutir uma série de métodos destinados à resolução de equações não lineares Zeros de funções reais - Introdução Necessidade de resolução de equações do tipo f(x) = 0 +FV -FV +FH-FH Em cada nó : FH = 0 FV = 0 FEstruturas (Lei de Kirchhoff) R E i v = g(i) + - E - Ri – g(i) = 0 Circuitos Zeros de funções reais - Introdução xÎ é um zero da função f(x) ou raiz da equação f(x) = 0 se f(x) = 0. Raízes podem ser números reais ou complexos. Trataremos somente de raízes reais de f(x). Abscissas dos pontos onde a curva intercepta o eixo x x2x2x1x1 f(x) x 18/02/2019 2 Zeros de funções reais - Introdução Para uma equação de segundo grau na forma: Determinação das raízes em função de a, b e c: Polinômios de grau mais elevado e funções com maior grau de complexidade Impossibilidade de determinação exata dos zeros Uso de soluções aproximadas 02 =++ cbxax a acbbx 2 42 --= Zeros de funções reais - Introdução Etapas para a determinação de raízes a partir de métodos numéricos FASE 1: Determinação de um intervalo (o menor possível) que contenha apenas uma raiz FASE 2: Melhoramento do valor da raiz aproximada (refinamento até que a raiz esteja dentro uma precisão ε prefixada) Zeros Reais de Funções Reais – Isolamento de Raízes Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Isolamento de raízes Realiza-se uma análise teórica e gráfica da função f(x). A precisão das análises é relevante para o sucesso da fase posterior. Teorema 1: Sendo f(x) contínua em um intervalo [a, b], se f(a)f(b) < 0 então existe pelo menos um ponto x = x entre a e b que é zero de f(x). Isolamento de raízes – Análise Gráfica xx11 xx22 f(x) xxx33 aa bb xx bb f(x) x aa aa xx11 f(x) xxx22 bb Isolamento de raízes – Tabelamento Exemplo: f(x) é contínua para I1 = [-5, -3] I2 = [0, 1] I3 = [2, 3] Cada um dos intervalos acima contém pelo menos um zero de f(x). Î x 39)( 3 +-= xxxf 18/02/2019 3 Isolamento de raízes – Tabelamento Exemplo: f(x) admite pelo menos um zero no intervalo [1,2] Mas esse zero é único? Análise do sinal de f’(x) f(x) admite um único zero em todo seu domínio de definição, localizado no intervalo [1,2] xexxf --= 5)( 0,05 2 1 )(' >>+= - xe x xf x Isolamento de raízes A partir do Teorema 1, se f’(x) existir e preservar o sinal em (a,b), então esse intervalo contém um único zero de f(x) Isolamento de raízes Se f(a)f(b) > 0, então se pode ter diversas situações no intervalo [a, b]. Isolamento de raízes A análise gráfica é fundamental para se obter boas aproximações para a raiz. Suficiente utilizar um dos seguintes passos: Esboçar o gráfico de f(x). Localizar as abscissas dos pontos onde a curva intercepta o eixo x. Obtenção da equação equivalente g(x) = h(x) a partir da equação f(x) = 0. Construção dos gráficos de g(x) e h(x) no mesmo sistema cartesiano e localização dos pontos x nos quais g(x) e h(x) se interceptam (f(x) = 0 g(x) = h(x)). Pode-se (deve-se) usar programas para traçar gráficos de funções dentro do intervalo de interesse. Isolamento de raízes O esboço do gráfico de uma função requer um estudo detalhado de seu comportamento, no qual devem ser considerados os itens abaixo: Domínio da função Pontos de descontinuidade Intervalos de crescimento e decrescimento Pontos de máximo e mínimo Concavidade Pontos de inflexão, etc Isolamento de raízes Exemplo: Solução utilizando o método 1: 39)( 3 +-= xxxf 30)(' 93)(' 39)( 2 3 == -= +-= xxf xxf xxxf )3,4(1 --Îx )1,0(2 Îx )3,2(3 Îx 33 -72 -7,3923 3 -51 30 11-1 13,3923- 3 3-3 -25-4 f(x)x x3x3 f(x) x-4 1-3 -2 -1 2 3 4 x2x2x1x1 18/02/2019 4 Isolamento de raízes Exemplo: Solução utilizando o método 2: Dada: Equação Equivalente: 039)( 3 =+-= xxxf 0393 =+- xx 3)( xxg = 39)( -= xxh )3,4(1 --Îx )1,0(2 Îx )3,2(3Îx x3x3 g(x) x-4 1-3 -2 -1 2 3 4x2x2 x1x1 h(x) y Isolamento de raízes Exemplo: Solução utilizando o método 2: Dada: Equação Equivalente: 05)( =-= - xexxf xex -= 5 xxg =)( xexh -= 5)( )2,1(Îx xx g(x) x1 2 3 4 h(x) y 5 6 Isolamento de raízes Exemplo: Solução utilizando o método 2: Dada: Equação Equivalente: 01)log()( =-= xxxf x x 1 )log( = )log()( xxg = x xh 1 )( = )3,2(Îx xx g(x) x1 2 3 4 h(x) y 5 6 Zeros Reais de Funções Reais – Refinamento de Raízes Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Refinamento de raízes Aplicação de métodos numéricos destinados ao refinamento de raízes: I. Método da Bisseção II. Método da Posição Falsa III. Método do Ponto Fixo IV. Método de Newton-Raphson V. Método da Secante Diferenciação dos métodos Modo de refinamento. Método Iterativo Caracterizado por uma série de instruções executáveis seqüencialmente, algumas das quais repetidas em ciclos (iterações). Refinamento de raízes Sequência de passos: 18/02/2019 5 Critérios de Parada Teste: xk suficientemente próximo da raiz exata? Como verificar tal questionamento? Interpretações para raiz aproximada x é raiz aproximada com precisão e se: ou Como proceder se não se conhece x ? ex <-x e<)(xf Critérios de Parada Reduz-se o intervalo que contém a raiz a cada iteração. Obtém um intervalo [a,b] tal que: . então Logo pode ser tomado como ï þ ï ý ü <- Î e x ab e ba, ex <-Î xbax ,, bax ,Î x Critérios de Parada Nem sempre é possível satisfazer ambos os critérios Critérios de Parada Métodos numéricos devem satisfazer a pelo menos um dos critérios Quando da utilização de programas computacionais, devemos utilizar: Teste de Parada Estipular o número máximo de iterações Prevenção de loops por: Erro no programa Escolha de método inadequado Zeros Reais de Funções Reais – Método da Bisseção Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Método da Bisseção Dada uma função f(x) contínua no intervalo [a,b] onde existe uma raiz única, é possível determinar tal raiz subdividindo sucessivas vezes o intervalo que a contém pelo ponto médio de a e b. Em outras palavras, o objetivo deste método é reduzir a amplitude do intervalo que contém a raiz até atingir precisão requerida, ou , usando para isto a sucessiva divisão de [a,b] ao meio e<- kk ab e<)(xf 18/02/2019 6 Método da Bisseção Definição do intervalo inicial Atribui-se [a,b] como intervalo inicial a0 = a b0 = b Condições de Aplicação f(a) x f(b) < 0 Sinal da derivada constante Método da Bisseção Definição de novos intervalos Calcula-se o ponto médio entre a e b, chamado de x0 Determina-se qual o subintervalo – [a , x0] ou [x0 , b] – contém a raiz Calcula-se o produto f(a) * f(x0) Verifica-se f(a) * f(x0) < 0 Se verdadeiro Logo a = a e b = x0 Caso contrario Logo a = x0 e b = b Repete-se o processo até que o valor de x atenda às condições de parada. ),( 0xaÎx ),( 0 bxÎx Método da Bisseção - Resumo ï þ ï ý ü ï î ï í ì > > < + = 0)( 0)( 0)( 2 0 0 0 00 0 xf bf af ba x ï î ï í ì = = Î Þ 01 01 00 ),( xb aa xax ï þ ï ý ü ï î ï í ì < > < + = 0)( 0)( 0)( 2 1 1 1 11 1 xf bf af ba x ï î ï í ì = = Î Þ 12 12 11 ),( bb xa bxx ï þ ï ý ü ï î ï í ì < > < + = 0)( 0)( 0)( 2 2 2 2 22 2 xf bf af ba x ï î ï í ì = = Î Þ 23 23 22 ),( bb xa bxx Método da Bisseção - Graficamente ba x0|| a1 x1 || a3 a2 || b1 || x2 || b3 x y b2= Método da Bisseção Exemplo: Utilizando o método de Equações Equivalentes para Isolamento de Raízes: Equação Equivalente: 01)log()( =-= xxxf x x 1 )log( = )log()( xxg = x xh 1 )( = )3,2(Îx h(x) y xx g(x) x1 2 3 4 5 6 Método da Bisseção Exemplo: 01)log()( =-= xxxf x0 = 2+3 2 = 2.5 f (2) = -0.3979 < 0 f (3)= 0.4314 > 0 f (2.5) = -5.1510-3 < 0 ì í ï î ï ü ý ï þ ï Þ x Î (2.5, 3) a1 = x0 = 2.5 b1 = b0 = 3 ì í ï î ï x1 = 2.5+3 2 = 2.75 f (2.5) < 0 f (3) > 0 f (2.75)= 0.2082 > 0 ì í ï î ï ü ý ï þ ï Þ x Î (2.5, 2.75) a2 = a1 = 2.5 b2 = x1 = 2.75 ì í ï î ï 18/02/2019 7 Método da Bisseção - Algoritmo k = 0; a0 = a; b0 = b; xk = (ak + bk)/2; while and if f(ak)f(xk) < 0 then /*raiz em [ak , xk] */ ak+1 = ak; bk+1 = xk; else /* raiz em [xk, bk] */ ak+1 = xk; bk+1 = bk ; end if xk+1 = (ak+1 + bk+1)/2; k = k +1; end while bk - ak >=e f (xk ) >= e Método da Bisseção - Algoritmo Ao final da execução do algoritmo, teremos um intervalo [ak, bk] que contém a raiz e uma aproximação para a raiz exata (tal que ou ) A convergência do método é intuitiva e<- kk ab x e<)(xf Método da Bisseção – Estimativa do número de iterações Após n iterações, a raiz estará contida no intervalo: Deve-se obter o valor de k, tal que , ou seja: 2 11 -- -=- kkkk ab ab e<- kk ab e< - k ab 2 00 )2log( )log()log( 00 e--> ab k k ab 2 00 -= e 002 abk -> )log()log()2log( 00 e--> abk Exemplo: Considerando um intervalo [2,3] e ε=10-2, calcular o numero mínimo de iterações para que tenhamos (Critério de Parada).e<- kk ab )2log( )log()log( 00 e--> ab k )2log( )10log()23log( 2--- >k 64,6 3010,0 2 )2log( )10log(2)1log( = + >k 7=k Método da Bisseção – Estimativa do número de iterações Método da Bisseção Exemplo: Utilizando o método de Equações Equivalentes para Isolamento de Raízes Equação Equivalente 039)( 3 =+-= xxxf 0393 =-= xx 3)( xxg = 39)( -= xxh )3,4(1 --Îx )1,0(2 Îx )3,2(3 Îx Método da Bisseção Exemplo: Cálculo da 1ª aproximação x0 = (a0+b0)/2 = (0+1)/2 = x0 = 0,5 f(x0) = 0,5 3 – 9x0,5 + 3 = -1,375 Teste de Parada |b-a| = |1| > 10-3 e |f(x0)| = |-1,375| = 1,375 > 10 -3 Escolha do Novo Intervalo f(a0) = 0 3 – 9x0 + 3 = 3, logo f(a0) > 0 f(b0) = 1 3 – 9x1 + 3 = -5, logo f(b0) < 0 f(x0) = 0,5 3 – 9x0,5 + 3 = -1,375, logo f(x0) < 0 logo: a1=a0=0 e b1=x0=0,5 039)( 3 =+-= xxxf 1,0=I 3103 -=e 18/02/2019 8 Método da Bisseção Exemplo: Então em 9 iterações . foi atendida, enquanto , não foie<- kk abe<)(xf 039)( 3 =+-= xxxf 1,0=I 3103 -=e 337890625.0=x Método da Bisseção Vantagens: Facilidade de implementação; Estabilidade e convergência para a solução procurada; Desempenho regular e previsível. Desvantagens Lentidão do processo de convergência (requer o cálculo de f(x) em um elevado número de iterações); Necessidade de conhecimento prévio da região na qual se encontra a raiz de interesse (nem sempre é possível); Complexidade da extensão do método para problemas multivariáveis. 258.24 10 3 7 00 =ÞÞ þ ý ü = =- - kk ab e Método da Bisseção – Exercício a) Analise a função abaixo. Identifique um intervalo onde existe raiz real. Explique porque essa raiz é única. Execute as primeiras 7 iterações do Método da Bisseção para a função , tal que b) Caso a condição de erro não tenha sido satisfeita, calcule quantas iterações ainda seriam necessárias. 1)( 3 --= xxxf 3102 -<e x1 2 3 4 y 50-1-2-3-4 1 2 3 4 -4 -3 -2 -1 Método da Bisseção – Exercício a) Execute as primeiras 5 iterações do Método da Bisseção para a função , tal que Para a iteração 5 temos: e 1)( 3 --= xxxf 3102 -<e Iter. a b f(a) f(b) x f(x ) 1 1,000000 2,000000 -1,000000 5,000000 1,500000 0,875000 2 1,000000 1,500000 -1,000000 0,875000 1,250000 -0,296875 3 1,250000 1,500000 -0,296875 0,875000 1,375000 0,224609 4 1,250000 1,375000 -0,296875 0,224609 1,312500 -0,051514 5 1,312500 1,375000 -0,051514 0,224609 1,343750 0,082611 31020,06253125,1375,1 -=-=- ab 31020,082611)( -=xf Método da Bisseção – Exercício b) Caso a condição de erro não tenha sido satisfeita, calcule quantas iterações ainda seriam necessárias. )2log( )log()log( 00 e--> ab k )2log( )102log()12log( 3--- >k )2log( )10log32(log)1log( -- >k 9658,8 30103,0 2,69897 30103,0 )330103,0(0 )2log( )10log32(log)1log( == -- = -- >k 9=k Zeros Reais de Funções Reais – Método da Posição Falsa Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira 18/02/2019 9 Método da Posição Falsa Método da Bisseção Calcula a média aritmética dos limites do intervalo que contém a raiz ([a, b]) Método da Posição Falsa Calcula a média ponderada dos limites do intervalo que contém a raiz ([a, b]) Método da Posição Falsa Dada a função e, sendo o intervalo inicial , temos que está mais próximo de zero que Logo é provável que a raiz esteja mais próxima de x = 0 que de x = 1 ( isso é sempre verdade quando f(x) é linear em ) Assim, ao invés de tomar a média aritmética, o método da posição falsa toma a média ponderada, com pesos de e 039)( 3 =+-= xxxf 1,0, =ba )0(305)1( ff =<<-= )0(f )1(f ba, )(af )(bf )()( )()( )()( )()( afbf abfbaf afbf afbbfa x - - = + + = Método da Posição Falsa - Graficamente Graficamente x é a interseção entre o eixo x e a reta que passa pelos pontos (a, f(a)) e (b, f(b)): Método da Posição Falsa - Graficamente x a = a0 xx f(x) b = b0x0 x0 = a0f(b0) - b0f(a0) f(b0) - f(a0) x a = a1 xx f(x) b1 = x1 x1 = a1f(b1) – b1f(a1) f(b1) - f(a1) x1 Método da Posição Falsa Definição do intervalo inicial Atribui-se [a,b] como intervalo inicial a0 = a b0 = b Condições de Aplicação f(a) x f(b) < 0 Sinal da derivada constante Método da Posição Falsa Definição dos Subintervalos Subdivide-se o intervalo pelo ponto de interseção da reta que liga f(a) a f(b) e o eixo das abscissas Verifica-se se, através do teste de parada, se x0 é uma aproximação da raiz da equação (x) pelo tamanho do intervalo [a, b] ou o valor f(x0) Se verdadeiro x0 é a raiz procurada Caso contrário define-se um novo intervalo 18/02/2019 10 Método da Posição Falsa Definição do novo intervalo Determina-se qual o subintervalo – [a , x0] ou [x0 , b] – contém a raiz Calcula-se o produto f(a) * f(x0) Verifica-se f(a) * f(x0) < 0 Se verdadeiro Logo a = a e b = x0 Caso contrario Logo a = x0 e b = b Repete-se o processo até que o valor de x atenda às condições de parada. ),( 0xaÎx ),( 0 bxÎx Método da Posição Falsa Exemplo: logo, existe ao menos 1 raiz no intervalo dado . Como têm o mesmo sinal, ]3,2[,1)log()( =-= Ixxxf þ ý ü >= <-= 04314,0)( 03979,0)( 0 0 bf af )()( )()( 0 afbf abfbaf x - - = )3979,0(4314,0 )3979,0(34314,02 -- -- = 8293,0 0565,2 = 4798,2= 00219,0)( 0 <-=xf )()( 00 xfeaf þ ý ü î í ì >== <== 0)(3 0)(4798,2 101 101 bfbb afxa Exemplo: Como , temos: Método da Posição Falsa ]3,2[,1)log()( =-= Ixxxf þ ý ü î í ì >= <== 0)(3 0)(4798,2 11 101 bfb afxa )0219,0(4314,0 )0219,0(34314,04798,2 -- -- = 0,4533 1354,1 = 5049,21 =x 00011,0)( 1 <-=xf þ ý ü î í ì >== <== 0)(3 0)(5049,2 112 112 bfbb afxa Método da Posição Falsa - Algoritmo k = 0; ak = a; bk = b; FAk = f(ak); GBk = f(bk); xk = (akGBk - bkFAk) / (GBk - FAk); while and if f(ak)f(xk) ≤ 0 then /* raiz em [ak , xk] */ bk = xk; else /* raiz em [xk, bk] */ ak = xk; end if FAk = f(ak); GBk = f(bk); xk = (akGBk - bkFAk) / (GBk - FAk); k = k +1; end while e>=- kk ab e>=)( kxf Método da Posição Falsa Exemplo: Cálculo da 1ª aproximação Teste de Parada |b-a| = |1| > 10-3 e |f(x0)| = |-0,322265625| > 10 -3 Escolha do Novo Intervalo f(a0) = 0 3 – 9x0 + 3 = 3, logo f(a0) > 0 f(b0) = 1 3 – 9x1 + 3 = -5, logo f(b0) < 0 f(x0) = 0,375 3 – 9x0,375 + 3 = -0,32..., logo f(x0) < 0 logo: a1=a0=0 e b1=x0=0,375 039)( 3 =+-= xxxf 1,0=I 3102 -=e )()( )()( 0 afbf abfbaf x - - = )3(5 )3(1)5(0 -- -- = 8 3 - - = Método da Posição Falsa Exemplo: Então em 3 iterações . foi atendida, enquanto , não foi No método da Bisseção, o valor foi encontrado depois de 9 iterações e<- kk abe<)(xf 039)( 3 =+-= xxxf 1,0=I 3103 -=e 337890625.0=x 337635046.0=x 18/02/2019 11 Método da Posição Falsa Vantagens: Estabilidade e convergência para a solução procurada; Desempenho regular e previsível; Cálculos mais simples que o método de Newton. Desvantagens: Lentidão do processo de convergência (requer o cálculo de f(x) em um elevado número de iterações); Necessidade de conhecimento prévio da região na qual se encontra a raiz de interesse (o que nem sempre é possível). Método da Posição Falsa– Exercício a) Analise a função abaixo. Identifique um intervalo onde existe raiz real. Execute as iterações do Método da Posição Falsa para a função , tal que 1)( 3 --= xxxf 3102 -<e x1 2 3 4 y 50-1-2-3-4 1 2 3 4 -4 -3 -2 -1 Recordando... O que é o Cálculo Numérico? O Cálculo Numérico corresponde a um conjunto de ferramentas ou métodos usados para se obter a solução de problemas matemáticos de forma aproximada, mas com um grau crescente de exatidão Qual o papel (ou função) do Cálculo Numérico na Engenharia? Solucionar problemas técnicos através de métodos numéricos, usando um modelo matemático Recordando... Em que consiste obter a raiz de uma função? Quais foram os métodos estudados até o momento para obtenção de raiz? Os métodos são capazes de obter todas as raízes da função em estudo? Existe alguma condição para a aplicação dos métodos estudados até o momento? Zeros Reais de Funções Reais – Método do Ponto Fixo Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Método do Ponto Fixo Dada uma função f(x) contínua no intervalo [a,b] onde existe uma raiz única, f(x) = 0, é possível transformar tal equação em uma equação equivalente x = φ(x) e, a partir de uma aproximação inicial x0, gerar uma sequência {xk} de aproximações para x pela relação xk+1 = φ(xk), uma vez que φ(x) é tal que f(x) = 0 se e somente se φ(x) = x. Transformamos o problema de encontrar zero de f(x) no problema de encontrar um ponto fixo de φ(x) A função φ(x) é chamada de função de iteração 18/02/2019 12 Método do Ponto Fixo Exemplo: Dada a função São funções de iteração possíveis: A forma geral das funções de iteração φ(x) é dada por com a condição de que A(x) 0 em x, ponto fixo de φ(x) 06)( 2 =-+= xxxf 1 6 )( 1 6 )( 6)( 6)( 4 3 2 2 1 + = -= -= -= x x x x xx xx )()()( xfxAxx += Método do Ponto Fixo A partir da definição da forma de φ(x), , podemos então mostrar que Existem infinitas equações de iteração φ(x) para a equação f(x) = 0 xxx == )(0)(f 0)( =Þ xx fquetalseja )()()( xxxx fA+= )0)(()( ==Þ xxx fporque xx = )(se xxxx =+Þ )()( fA 0)()( =Þ xx fA )0)((0)( =Þ xx Aporquef )()()( xfxAxx += Método do Ponto Fixo - Graficamente Uma raiz da equação φ(x)=x é a abscissa do ponto de interseção da reta y=x com a curva y=φ(x) Método do Ponto Fixo - Graficamente Todavia, para algumas escolhas de φ(x) o Método do Ponto Fixo pode divergir do valor x procurado Método do Ponto Fixo Exemplo: Dada a equação : As raízes são x1 = -3 e x2 = 2 (Não há necessidade de uso de métodos numéricos para o calculo) Objetivo: Mostrar a convergência ou divergência do processo iterativo Seja a raiz x2 = 2 e φ1 (x) = 6 - x 2,Tomando x0= 1,5 e φ (x) = φ1 (x) Seja a raiz x2 = 2 e ,Tomando x0= 1,5 e φ (x) = φ2 (x) 06)( 2 =-+= xxxf xx -= 6)(2 Método do Ponto Fixo Exemplo: Dada a equação , com raiz x2 = 2 , φ1 (x) = 6 - x 2 e x0 = 1,5 x1 = φ(x0) = 6 – 1,5 2 = 3,75 x2 = φ(x1) = 6 – 3,75 2 = -8,0625 x3 = φ(x2) = 6 – (-8,0625) 2 = -59,003906 x4 = φ(x3) = 6 – (-59,003906) 2 = -3475,4609 Conclui-se que {xk} não convergirá para x2 = 2 06)( 2 =-+= xxxf 18/02/2019 13 Método do Ponto Fixo 06)( 2 =-+= xxxf Exemplo: Dada a equação , com raiz x2 = 2 , φ1 (x) = 6 - x 2 e x0 = 1,5 Análise Gráfica: y xx2x2 x1 φ(x) xx00 y = x x2 x 1 x 1 {xk} x Método do Ponto Fixo Exemplo: Dada a equação , com raiz x2 = 2 , e x0 = 1,5 x1 = φ(x0) = x2 = φ(x1) = x3 = φ(x2) = x4 = φ(x3) = x5 = φ(x4) = Conclui-se que {xk} tende a convergir para x2 = 2 06)( 2 =-+= xxxf xx -= 6)(2 121320343,25,16 =- 969436380,1121320343,26 =- 007626364,2969436380,16 =- 998092499,1007626364,26 =- 000476818,2998092499,16 =- Método do Ponto Fixo Exemplo: Dada a equação , com raiz x2 = 2 , e x0 = 1,5 Análise Gráfica: 06)( 2 =-+= xxxf xx -= 6)(2 {xk} x2 quando k inf φ(x) x y y = x x 2 x 2 x1 x0 x2 Método do Ponto Fixo Teorema 2: Sendo x uma raiz de f(x) = 0, isolada em um intervalo I centrado em x e φ(x) uma função de iteração para f(x) = 0. Se 1. φ(x) e φ’(x) são contínuas em I 2. |φ’(x)| < 1, x Î I e 3. x0 Î I então a sequencia {xk} gerada pelo processo iterativo xk+1 = φ(xk) convergirá para x . Além disso quanto menor for o valor de |φ’(x)|, mais rápido o Método do Ponto Fixo convergirá. Método do Ponto Fixo Resgatando os exemplos anteriores, para a função temos que: φ1(x) ( ) geração de uma sequencia divergente de x2 = 2 φ2(x) ( ) geração de uma sequencia convergente para x2 = 2 06)( 2 =-+= xxxf 2 1 6)( xx -= xx -= 6)(2 Método do Ponto Fixo Avaliando a divergência de φ1(x) φ1(x) = 6 - x 2 e φ’1(x) = - 2x contínuas em I |φ’1 (x)| < 1 |-2x| < 1 -½ < x < ½ Não existe um intervalo I centrado em x2=2, tal que |φ’(x)| < 1, x Î I φ1 (x) não satisfaz a condição 2 do Teorema 2 com relação a x2=2. 18/02/2019 14 Método do Ponto Fixo Avaliando a convergência de φ2(x) e φ2 (x) é contínua em S = {x Î R | x 6} φ’2 (x) é contínua em S’ = {x Î R | x < 6} . É possível obter um intervalo I centrado em x2=2, tal que todas as condições do Teorema 2 sejam satisfeitas. xx -= 6)(2 )62/1()('2 xx --= 75,5162/11)('2 <<-< xxx Método do Ponto Fixo - Algoritmo Critérios de Parada |f(xk)| e |xk – xk-1| e k = 0; Xk+1 = φ(xk); while and k = k +1; xk = xk+1; xk+1 = φ(xk); end while e>=-+ kk xx 1 e>=+ )( 1kxf Método do Ponto Fixo – Verificando a Convergência Exemplo: Dada a função , cujas raízes são 2 e -3, vamos avaliar a convergência da função equivalente , dados x1 = -3 e x0= -2,5 06)( 2 =-+= xxxf 1 6 )(3 -= x x 0,,0 6 )(' 2 Î< - = xx x x 0,, 66 )(' 22 Î= - = xx xx x 6661 6 1)(' 2 2 >-<><< xouxx x x 0,,1 6 )( Î-= xx x x Método do Ponto Fixo – Verificando a Convergência Exemplo: Dada a função , cujas raízes são 2 e -3, vamos avaliar a convergência da função equivalente , dados x1 = -3 e x0= -2,5 Como o objetivo é obter a raiz negativa, temos: Podemos então definir o intervalo que o processo convergirá visto que o intervalo está centrado na raiz x = -3 )6;(:,,1)(' 111 --=Î< IseráIxxquetalI )4497897,26( -=- )5.2,5.3( --=I 1II 06)( 2 =-+= xxxf 1 6 )(3 -= x x Método do Ponto Fixo – Verificando a Convergência Exemplo: Dada a função , cujas raízes são 2 e -3, vamos avaliar a convergência da função equivalente , dados x1 = -3 e x0= -2,5 Tomando x0= -2,5, temos: Quando não se conhece a raiz, escolhe-se o intervalo I aproximadamente centrado em x Quanto mais preciso isolamento de x, maior exatidão na escolha de I 892617,2 170213,3 764706,2 5,2 4 3 2 1 -= -= -= -= x x x x 06)( 2 =-+= xxxf 1 6 )(3 -= x x Método do Ponto Fixo Exemplo: Dados: , calcule a raiz de f(x) utilizando o MPF: Assim, e Importante lembrar: Iteramos de modo que , todavia avaliamos, a cada iteração se Desafio: Provar que satisfaz a condição 2 do Teorema 2 no intervalo (0, 1) ; 3 1 9 )(;039)( 3 3 +==+-= x xxxxf )1,0(;105;5,0 20 Î== - xex 3376233,0=x 31012,0)( --=xf )(1 kk xx =+ e<)( kxf )(x 18/02/2019 15 Método do Ponto Fixo Vantagens Rapidez processo de convergência; Desempenho regular e previsível. Desvantagens Um inconveniente é a necessidade da obtenção de uma função de iteração φ(x); Difícil sua implementação. Método do Ponto Fixo – Exercício Resolvido 1) Tente encontrar a raiz da função utilizando a função de iteração e , sendo . Analise sua resposta. 1)( 3 --= xxxf 3102 -<e x1 2 3 4 y 50-1-2-3-4 1 2 3 4 -4 -3 -2 -1 2 11 )( xx x += 10 =x Método do Ponto Fixo – Exercício Resolvido 1) Tente encontrar a raiz da função utilizando a função de iteração e , sendo x1 = φ(x0) = x2 = φ(x1) = x3 = φ(x2) = x4 = φ(x3) = x5 = φ(x4) = 1)( 3 --= xxxf 3102 -<e 2 11 )( xx x += 10 =x 2 1 1 1 1 2 =+ 75,0 2 1 2 1 2 =+ ...1111,3 75,0 1 75,0 1 2 =+ ...4247,0 1111,3 1 1111,3 1 2 =+ ...8973,7 4247,0 1 4247,0 1 2 =+ Método do Ponto Fixo – Exercício Resolvido 1) Tente encontrar a raiz da função utilizando a função de iteração e , sendo x6 = φ(x5) = x7 = φ(x6) = Conclui-se que {xk} tende a divergir da raiz da equação f(x). 1)( 3 --= xxxf 3102 -<e 2 11 )( xx x += 10 =x ...1427,0 8973,7 1 8973,7 1 2 =+ ...1461,56 1427,0 1 1427,0 1 2 =+ Método do Ponto Fixo – Exercício Resolvido 1) Tente encontrar a raiz da função utilizando a função de iteração e , sendo Justificando a resposta: Como a condição deve ser satisfeita, onde I é o intervalo centrado em x , é fácil perceber que isso não acontece, uma vez que 1)( 3 --= xxxf 3102 -<e 2 11 )( xx x += 10 =x 0, 11 )( 2 Î+= xx xx x 0, 21 )(' 32 Î- - = xx xx x 1 2 1 2 1 21 1)(' 33332 < -- <- - <- - < x x xx x xx x Ixx Î<1)(' 03)1(')(' 00 >==Î xeIx Método do Ponto Fixo – Exercício 1) Tente encontrar a raiz da função utilizando a função de iteração e , sendo .Justifique sua resposta. 1)( 3 --= xxxf 3102 -<e x1 2 3 4 y 50-1-2-3-4 1 2 3 4 -4 -3 -2 -1 13 1 )( 2 3 - -- -= x xx xx 10 =x 18/02/2019 16 Zeros Reais de Funções Reais – Método de Newton Raphson Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Método de Newton-Raphson Método do Ponto Fixo (MPF) Uma das condições de convergência é que |φ’(x)| M < 1, x Î I , onde I é um intervalo centrado na raiz A convergência será tanto mais rápida quanto menor for |φ’(x)| O método de Newton busca garantir e acelerar a convergência do MPF Escolha de φ(x), tal que φ’(x) = 0, como função de iteração Método de Newton-Raphson Dada a equação f(x) = 0 e partindo da forma geral para φ(x) φ(x) = x + A(x)f(x) Busca-se obter a função A(x) tal que φ’(x) = 0 φ(x) = x + A(x)f(x) Þ φ’(x) = 1 + A’(x)f(x) + A(x)f’(x) Þ φ’(x) = 1 + A’(x)f(x) + A(x)f’(x) Þ φ’(x) = 1 + A(x)f’(x) Método de Newton-Raphson Assim donde se toma Como φ(x) = x + A(x)f(x) Logo: )( )(' 1 )( xf xf xx - += -= )(' )( )( xf xf xx 0)(' =x 0)(')(1 =+ xx fA )(' 1 )( x x f A - = )(' 1 )( xf xA - = Método de Newton-Raphson Então, dada f(x), a função de iteração φ(x) = x - f(x)/f’(x) será tal que φ’(x) = 0, posto que e, como f(x) = 0, φ’(x) = 0 ( desde que f’(x) 0 ) - -= 2 2 )]('[ )('')()]('[ 1)(' xf xfxfxf x 2 2 2 2 )]('[ )('')()]('[ )]('[ )]('[ )(' xf xfxfxf xf xf x - -= 2)]('[ )('')( )(' xf xfxf x = Método de Newton-Raphson Deste modo, escolhido x0, a sequência {xk} será determinada por onde k = 0, 1, 2, ... )(' )( 1 k k kk xf xf xx -=+ 18/02/2019 17 Método de Newton-Raphson - Convergência Teorema 3: Sendo f(x), f’(x) e f”(x) contínuas em um intervalo I que contém uma raiz x = x de f(x) = 0 e supondo f’(x) 0, existirá um intervalo Ī I contendo a raiz x, tal que se x0 Î Ī, a seqüência {xk} gerada pela fórmula recursiva convergirá para a raiz. )(' )( 1 k k kk xf xf xx -=+ Método de Newton-Raphson – Graficamente Dado o ponto ( xk , f(xk) ) Traçamos a reta Lk(x) tangente à curva neste ponto: Lk(x) = f(xk) + f’(xk)(x-xk) Determinanos o zero de Lk(x), que é um modelo linear que aproxima f(x) em uma vizinhança xk Faz-se xk +1 = x 0)( =xLk )(' )( k k k xf xf xx -= Método de Newton-Raphson – Graficamente Análise Gráfica Repete-se o processo até que o valor de x atenda às condições de parada. x xxxx f(x) x1xx00 x2 x3 1a iteração 2a iteração 3a iteração 4a iteração Método de Newton-Raphson - Algoritmo Teste de parada: |f(xk)| ε |xk – xk-1| ε Algoritmo: x0 := x; k := 0; while |f(xk)| > ε and |xk – xk-1| > ε xk+1 := xk – f(xk)/f’(xk) k := k +1; end while Método de Newton-Raphson Exemplo: Dado f(x) = x2 + x – 6 , x2 = 2 e x0 = 1,5 Fórmula recursiva: 12 6 )(' )( )( 2 + -+ -=-= x xx x xf xf xx 0625,2 15,12 65,15,1 5,1)( 2 01 = + -+ -== xx 000762195,2 10625,22 60625,20625,2 0625,2)( 2 12 = + -+ -== xx 000000116,2)( 23 == xx Método de Newton-Raphson Exemplo: Dado f(x) = x2 + x – 6 , x2 = 2 e x0 = 1,51,5 Comentários: A parada poderá ocorrer na 3a iteração (x = 2,000000116), caso a precisão do cálculo com 6 casas decimais seja satisfatória para o contexto do trabalho Observe que, no Método do Ponto Fixo, com o valor x = 2,000476818 foi encontrado somente na 5a iteração xx -= 6)( 18/02/2019 18 Método de Newton-Raphson Exemplo: Considere a função f(x) = x3 - x - 1, e ε = 0,002 cujos zeros encontram-se nos intervalos: x1 Î I1 = (-1, 0), x2 Î I2 = (1, 2) Seja: x0 = 1 )(' )( 1 k k kk xf xf xx -=+ 13 1 )( 2 3 - -- -= x xx xx Método de Newton-Raphson Exemplo: Considere a função f(x) = x3 - x - 1, e ε = 0,002 cujos zeros encontram-se nos intervalos: x1 Î I1 = (-1, 0), x2 Î I2 = (1, 2) Cálculo da 1ª aproximação φ(x0) = 1 – [ (1)³ – 1 – 1 ] = 1,5 = x1 [ 3x(1)² – 1 ] Teste de Parada |f(x1)| = | (1,5)³ – 1,5 – 1 | = 0,875 > e |x1-x0| =| 1,5 - 1 | = 0,5 > e Exemplo: Considere a função f(x) = x3 - x - 1, e ε = 0,002 cujos zeros encontram-se nos intervalos: x1 Î I1 = (-1, 0), x2 Î I2 = (1, 2) Cálculo da 2ª aproximação φ(x1) = 1,5 – [ (1,5)³ – 1,5 – 1 ] = 1,3478261 = x2 [ 3x(1,5)² – 1 ] Teste de Parada |f(x2)| = | 0,100682 | = 0,100682 > e |x2-x1| =| 1,3478261 - 1,5 | = 0,1521739 > e Método de Newton-Raphson Cálculo da 3ª aproximação φ(x2) = 1,3478261 - [ (1,3478261)³ - 1,3478261 - 1 ] [ 3x(1,3478261)² - 1 ] φ(x2) = 1,3252004 = x3 Teste de Parada |f(x3)| =| 0,0020584 | = 0,0020584 > e |x3-x2| =| 1,3252004 – 1,3478261 | = 0,0226257 > e Exemplo: Considere a função f(x) = x3 - x - 1, e ε = 0,002 cujos zeros encontram-se nos intervalos: x1 Î I1 = (-1, 0), x2 Î I2 = (1, 2) Método de Newton-Raphson Exemplo: Considere a função f(x) = x3 - x - 1, e ε = 0,002 cujos zeros encontram-se nos intervalos: x1 Î I1 = (-1, 0), x2 Î I2 = (1, 2) A sequência {xk} gerada pelo método de Newton será: Método de Newton-Raphson Iteração x |xk-xk-1| F(x) 1 1,5 0,5 0,875 2 1,3478261 0,1521739 0,1006822 3 1,3252004 0,0226257 0,0020584 4 1,3247182 0,0004822 1,0352x10-6 e = 0,002 Comprovando o impacto de uma boa escolha de x0 Exemplo: Considere a função f(x) = x3 – 9x + 3, que possui três zeros: x1 Î I1 = (-4, -3), x2 Î I2 = (0, 1) e x3 Î I3 = (2, 3). Seja x0 = 1,5: Método de Newton-Raphson 18/02/2019 19 Comprovando o impacto de uma boa escolha de x0 Exemplo: Considere a função f(x) = x3 – 9x + 3, que possui três zeros: x1 Î I1 = (-4, -3), x2 Î I2 = (0, 1) e x3 Î I3 = (2, 3). Seja x0 = 1,5: No início há um divergência da região onde estão as raízes, mas depois de x7 os valores se aproximam cada vez mais de x3 Causa: x0 (1,5) é próximo de , que é raiz de f´(x) Da mesma forma, x1 (-1,6666667) está próximo de , outra raiz de f’(x) Método de Newton-Raphson 3 3- Vantagens: Rapidez processo de convergência Desempenho elevado Desvantagens: Necessidade da obtenção de f’(x) , o que pode ser impossível em determinados casos O cálculo do valor numérico de f’(x) a cada iteração Método de Newton-Raphson Exercício Encontre uma raiz da função f(x) = x3-2x utilizando o método Newton-Raphson. Considere o gráfico a seguir: Zeros Reais de Funções Reais – Método da Secante Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Método da Secante Método de Newton-Raphson Um grande inconveniente é a necessidade da obtenção de f’(x) e o cálculo de seu valor numérico a cada iteração Forma de desvio do inconveniente Substituição da derivada f’(xk) pelo quociente das diferenças 1 1)()()(' - - - - kk kk k xx xfxf xf Método da Secante A função de iteração será: 1 1)()( )( )( - - - - -= kk kk k k xx xfxf xf xx 1 1)()( )( )( - - - - -= kk kk k k xx xfxf xf xx )()( )()( )()( )()( )( 1 1 1 1 - - - - - - - - - = kk kkkk kk kkkk xfxf xfxxfx xfxf xfxxfx x )()( )()( )( 1 11 - -- - - = kk kkkk xfxf xfxxfx x 18/02/2019 20 Método da Secante - Geometricamente A partir de duas aproximações xk-1 e xk obtém-se o ponto xk+1 como sendo a abscissa do ponto de intersecção do eixo x e da reta que passa pelos pontos ( xk-1 , f(xk-1) ) e ( xk , f(xk) ) (secante à curva da função) x 1a iteração 2a iteração 3a iteração 4a iteração xx f(x) x1xx00 x2 x3 x4 x5 Repete-se o processo até que o valor de x atenda às condições de parada. Método da Secante - Convergência Como o Método da Secante é uma aproximação do método de Newton, as condições de convergência são praticamente as mesmas, ou seja basta que o Teorema 3 seja satisfeito Todavia, o Método da Secante pode divergir para o seguinte caso )()( 1- kk xfxf )()( )()( )( 1 11 - -- - - = kk kkkk xfxf xfxxfx x Método da Secante - Algoritmo Testes de Parada |f(xk)| ε |xk – xk-1| ε Algoritmo x0 := x; x1 := x1; k := 1; x2 := (x0*f(x1) – x1*f(x0)) / (f(x1) - f(x0)); while |f(xk+1)| > ε and |xk+1 – xk| > ε xk+1 := (xk-1*f(xk) - xk*f(xk-1)) / (f(xk) - f(xk-1)); k := k +1; end while Método da Secante Exemplo: Considere-se a função f(x) = x3 - x - 1, e = 0,003, x0 = 1,5 e x1 = 1,7: )()( )()( )( 1 11 - -- - - = kk kkkk xfxf xfxxfx x Método da Secante Exemplo: Considere-se a função f(x) = x3 - x - 1, e = 0,003, x0 = 1,5 e x1 = 1,7: Cálculo da 1ª aproximação x0 = 1,5 e x1 = 1,7 f(x0) = 0,875 > 0 f(x1) = 2,213 > 0 x2 = [1,5 x (2,213) – 1,7 x (0,875)] = 1,36921 [2,213 – (0,875)] Teste de Parada |f(x2)| = | (1,36921)³ – 1,36921 – 1 | = 0,19769 > e |x2 - x1| =|1,36921 – 1,7| = 0,33079 > e Novo Intervalo: x1 = 1,7 e x2 = 1,36921 Método da Secante Exemplo: Considere-se a função f(x) = x3 - x - 1, e = 0,003, x0 = 1,5 e x1 = 1,7: Cálculo da 2ª aproximação x1 = 1,7 e x2 = 1,36921 f(x1) = 2,213 > 0 f(x2) = 0,19769 > 0 x3 = [1,7 x (0,19769) - 1,36921x (2,213)] = 1,33676 [0,19769 - 2,213] Teste de Parada |f(x3)| = |0,05193| = 0,05193 > e |x3 - x2| =|1,33676 – 1,36921| = 0,03245 > e Novo Intervalo: x2 = 1,36921 e x3 = 1,33676 18/02/2019 21 Método da Secante Exemplo: Considere-se a função f(x) = x3 - x - 1, e = 0,003, x0 = 1,5 e x1 = 1,7: Cálculo da 3ª aproximação x2 = 1,36921 e x3 = 1,33676 f(x2) = 0,19769 > 0 f(x3) = 0,05193 > 0 x4 = [1,36921 x (0,05193) - 1,33676 x (0,19769)] = [(0,05193) - 0,19769] x4 = 1,3252 Teste de Parada |f(x4)| = |0,00206| = 0,00206 < e cond. atendida |x4 – x3| =|1,3252 – 1,33676| = 0,01156 > e Método da Secante Exemplo: Considere-se a função f(x) = x2 + x – 6 = 0, x0 = 1,5 e x1 = 1,7: Solução: )()( )()( 01 0110 2 xfxf xfxxfx x - - = 25,241,1 25,27,1)41,1(5,1 +- -- = 2,035712 =x )()( )()( 12 1221 3 xfxf xfxxfx x - - = 1,99774= )()( )()( 23 2332 4 xfxf xfxxfx x - - = 1,99999= Método da Secante Exemplo: Considere-se a função f(x) = x2 + x – 6 = 0, x0 = 1,5 e x1 = 1,7: Comentários: A parada poderá ocorrer na 3a iteração (x = 1,99999 ), caso a precisão do cálculo com 5 casas decimais seja satisfatória para o contexto do trabalho No método de Newton-Raphson o valor x = 2,000000116, foi encontrado também na 3a iteração Método da Secante Vantagens Rapidez processo de convergência Cálculos mais convenientes que do método de Newton Desempenho elevado Desvantagens Se o cálculo f’(x) não for difícil, então o método logo será substituído pelo de Newton-Raphson Se o gráfico da função for paralela a um dos eixos e/ou tangencia o eixo das abscissas em um ou mais pontos, logo não se deve usar o Método da Secante Zeros Reais de Funções Reais – Comparação entre os métodos Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Comparação entre os métodos Critérios de Comparação Garantias de Convergência Rapidez de Convergência Esforço Computacional 18/02/2019 22 Comparação entre os métodos Garantias de Convergência Bissecção e Posição Falsa Convergência garantida, desde que a função seja contínua num intervalo [a,b] , tal que f(a)f(b) < 0 Ponto Fixo, Newton-Raphson e Secante Condições mais restritivas de convergência Se as condições de convergência forem satisfeitas, os dois últimos métodos são mais rápidos do que os demais estudados Comparação entre os métodos Rapidez de Convergência Número de Iterações Medida usualmente adotada para a determinação da rapidez de convergência de um método Não deve ser uma medida conclusiva sobre o tempo de execução do programa Tempo gasto na execução de uma iteração Variável de método para método Comparação entre os métodos Esforço Computacional Indicadores Número de operações efetuadas a cada iteração Complexidade das operações Número de decisões lógicas Número de avaliações de função a cada iteração Número total de iterações Comparação entre os métodos Esforço Computacional Conclusões gerais sobre a eficiência computacional de um método. Bissecção Cálculos mais simples por iteração Newton Cálculos mais elaborados Número de iterações da Bissecção é, na grande maioria das vezes, muito maior do que o número de iterações efetuadas por Newton Comparação entre os métodos Condições a Serem Satisfeitas pelo Método Ideal Convergência assegurada Ordem de convergência alta Cálculos por iteração simples Escolha do melhor método: Newton-Raphson Caso seja fácil a verificação das condições de convergência e o cálculo de f’(x) Secante Caso seja trabalhoso obter e/ou avaliar f’(x) , uma vez que não é necessária a obtenção de f’(x) Comparação entre os métodos Critério de Parada Detalhe importante na escolha do método: Se o objetivo for a redução do intervalo que contém a raiz Bissecção (não usar o Método da Posição Falsa) Se a escolha parte de um valor inicial para a raiz Newton-Raphson ou da Secante (pois trabalham com aproximações xk para a raiz exata) 18/02/2019 23 Comparação entre os métodos Observações importantes: Situações nas quais se deve evitar o uso do Método de Newton-Raphson e da Secante: Tendência da curva ao paralelismo a qualquer um dos eixos Tendência da função à tangência ao eixo das abscissas em um ou mais pontos. Comparação entre os métodos Conclusão: Escolha do método Diretamente relacionada com a equação cuja solução é desejada Comportamento da função na região da raiz exata Dificuldades com o cálculo de f´(x) Critério de parada, etc. Comparação entre os métodos Exemplo: f(x) = x3 – x – 1 x1 2 3 4 y 50-1-2-3-4 1 2 3 4 -4 -3 -2 -1 x Î [1, 2 ], e = 10-6 Comparação entre os métodos Exemplo: f(x) = x3 – x – 1 Observações: Melhor desempenho: Método do Ponto Fixo Métodos de Newton e Secante: muitas iterações pois houve divergência no inicio da execução ( denominador 0 ) Comparação entre os métodos Exemplo: f(x) = 4 sen(x) – e2 Observações: Melhor desempenho: Método de Newton, devido à boa escolha de x0 x Î [0, 1 ], e = 10-5 04_Sistemas_Lineares_Folheto.pdf 18/02/2019 1 Sistemas Lineares Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Programa 1. Introdução 2. Métodos Diretos a) Eliminação de Gauss b) Decomposição LU 3. Métodos Iterativos a) Gauss-Jacobi b) Gauss-Siedel Sistemas Lineares Introdução Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Introdução A resolução de sistemas lineares é um problema que surge nas mais diversas áreas Ex: Cálculos de estruturas, Redes de transporte, Redes de comunicação, etc Introdução Exemplo: Calcular tensões dos nós do circuito elétrico: Solução: Temos que a corrente entre 2 pontos é dada por: . Pela lei de Kirchoff a soma das correntes que chega a um nó é igual a soma das correntes que saem dele. Assim: 1V R VV I ba - = Introdução Exemplo: Calcular tensões dos nós do circuito elétrico: Nó 1: Nó 3: Nó 4: 1V 1 0 221 1141312 -= - + - + - VVVVVVV 31 2312 VVVV -= - 3 127 213 3134323 VVVVVVV -= - + - + - 21 1443 VVVV -= - Nó 2: 18/02/2019 2 Introdução Exemplo: Calcular tensões dos nós do circuito elétrico: Simplificando as equações: Nó 1: Nó 2: Nó 3: Nó 4: 1 0 221 1141312 -= - + - + - VVVVVVV 31 2312 VVVV -= - 3 127 213 3134323 VVVVVVV -= - + - + - 21 1443 VVVV -= - 026 4321 =+++- VVVV 043 321 =-+- VVV 25461323 4321 =-+-- VVVV 032 431 =-+ VVV Introdução Exemplo: Calcular tensões dos nós do circuito elétrico: Montando o sistema: Nosso problema agora se resume em encontrar os valores de V1, V2, V3 e V4 que solucionem o sistema linear acima. ï ï î ï ï í ì =-++ =-+-- =+-+- =+++- 03201 25461323 00143 01126 4321 4321 4321 4321 VVVV VVVV VVVV VVVV Introdução Um sistema linear com m equações e n variáveis tem a seguinte forma geral: onde: aij coeficientes 1 ≤ i ≤ m, 1 ≤ j ≤ n xj incógnitas j = 1,...,n bi termos independentes i = 1,...,m mnmnmm nn nn bxaxaxa bxaxaxa bxaxaxa =+++ =+++ =+++ ... ... ... 2211 22222121 11212111 Introdução Exemplo: onde: 2, 4, -5, 4, 1, -5, 2, 4 e 5 coeficientes x1, x2 e x3 incógnitas 5, 2 e -1 termos independentes ï î ï í ì -=++ =-+ =-+ 1542 2514 5542 321 321 321 xxx xxx xxx Introdução Forma Matricial: Ax = b na qual: ù é = mnmmm n n aaaa aaa aaa A 321 22221 11211 ù é = mb b b b 2 1 ù é = nx x x x 2 1 Introdução Exemplo: Forma Geral: Forma Matricial: ï î ï í ì -=++ =-+ =-+ 1542 2514 5542 321 321 321 xxx xxx xxx ù é - = ù é ´ ù é - - 1 2 5 542 514 542 3 2 1 x x x 18/02/2019 3 Introdução Relembrando... Multiplicação de Matrizes O produto de uma matriz A de dimensão n x m por um escalar k resulta em uma matriz B = kA de mesma dimensão n x m, tal que Ex: e .,, jikab ijij = ù é = 654 321 A ù é == 12108 642 2AB Introdução Relembrando... Multiplicação de Matrizes O produto de uma matriz A (n x m) por um vetor v (m x 1) resulta em um vetor x (n x 1) de forma que Ex: ....,,2,1, 1 nivax m j jiji == = ù é ==® ù é = ù é = 17 11 5 2 1 , 65 43 21 AvxvA Introdução Relembrando... Multiplicação de Matrizes O produto de uma matriz A (n x p) por uma matriz B (p x m) é uma matriz C = AB (n x m) tal que o elemento cij é obtido pela soma dos produtos da linha i de A pelos correspondentes elementos da coluna j de B. Logo, para a multiplicação de duas matrizes, o número de colunas da primeira tem que ser igual ao número de linhas da segunda Ex: ....,,2,1...,,2,1, 1 mjenivac p k kjikij === = 22 23 32 4841 126 53 04 61 , 653 012 x x x ABCBA ù é ==® ù é = ù é = Introdução – Classificação de Sistemas Classificação dos sistemas Solução Única det (A) ≠ 0 Infinitas Soluções ou Sem Solução det (A) = 0 Introdução – Classificação de Sistemas Relembrando…. Conceito de Determinante Uma matriz quadrada (n x n) A, chamada matriz de ordem n, tem um número associado, conhecido por determinante, cujo valor pode ser obtido pela fórmula de recorrência onde Mij é a matriz de ordem n-1 resultante da remoção da linha i e coluna j de A e sendo o determinante de uma matriz (1 x 1) igual a esse único elemento. Logo: )det()1(...)det()det()det( 11 1 12121111 nn n MaMaMaA +-++-= 1111 )det(][ aAaA =®= 12212211 2221 1211 )det( aaaaA aa aa A -=® ù é = Introdução – Classificação de Sistemas Relembrando…. Conceito de Determinante Matriz A com det(A) = 0 Matriz Singular Matriz A com det(A) ≠ 0 Matriz Não Singular ® ù é = 333231 232221 131211 aaa aaa aaa A )()()()det( 223132211323313321122332332211 aaaaaaaaaaaaaaaA -+---= 18/02/2019 4 Introdução – Classificação de Sistemas Solução Única Exemplo: det (A) = -6 -1 = -7 Infinitas Soluções Ex: det (A) = 4 - 4 = 0 î í ì -=- =+ 23 32 21 21 xx xx ù é =® 1 1 x î í ì =+ =+ 624 32 21 21 xx xx ù é - =® 23 x Introdução – Classificação de Sistemas Sem Solução Ex: det (A) = 4 - 4 = 0 2x1 + x2 =1 4x1 + 2x2 = 6 ì í î Introdução – Sistemas Triangulares Possibilidade de resolução da forma Direta Sistema Triangular Inferior Sistema Triangular Superior Introdução – Sistemas Triangulares Sistema Triangular Inferior A solução é calculada pelas substituições sucessivas , , ù é = ù é ´ ù é nnnnnnn b b b b x x x x aaaa aaa aa a 3 2 1 3 2 1 321 333231 2221 11 0 00 000 11 1 11111 a b xbxa =®= 22 1212 22222121 a xab xbxaxa - =®=+ 23 2321313 33333232131 a xaxab xbxaxaxa -- =®=++ Introdução – Sistemas Triangulares Sistema Triangular Inferior ù é = ù é ´ ù é nnnnnnn b b b b x x x x aaaa aaa aa a 3 2 1 3 2 1 321 333231 2221 11 0 00 000 nnnnnnnnnn bxaxaxaxa =++++ -- 11,2211 ... nn nnnnnn n a xaxaxab x 11,2211 ... ------=® Introdução – Sistemas Triangulares Sistema Triangular Inferior As substituições sucessivas podem ser representadas por: ù é = ù é ´ ù é nnnnnnn b b b b x x x x aaaa aaa aa a 3 2 1 3 2 1 321 333231 2221 11 0 00 000 ....,,2,1, 1 1 ni a xab x ii i j jiji i = - = - = 18/02/2019 5 Introdução – Sistemas Triangulares Sistema Triangular Inferior Exemplo: Calcular a solução do sistema triangular inferior usando as substituições sucessivas: , ù é = ù é ù é -- - 6 48 1 4 9341 0861 0053 0002 4 3 2 1 x x x x 2 2 4 42 11 ==®= xx 15 231 153 221 -= ´- =®=+ xxx 5 8 )1(6248 4886 3321 = -+- =®=+- xxxx Introdução – Sistemas Triangulares Sistema Triangular Inferior Exemplo: Calcular a solução do sistema triangular inferior usando as substituições sucessivas: Logo, o vetor solução é dado por: ù é = ù é ù é -- - 6 48 1 4 9341 0861 0053 0002 4 3 2 1 x x x x 3 9 )5(3)1(4)2(6 6934 44321 = +--+ =®=+-+- xxxxx ù é - = 3 5 1 2 x Introdução – Sistemas Triangulares Sistema Triangular Superior A solução é calculada pelas substituições retroativas: , ù é = ù é ´ ù é nnnn n n n d d d d x x x x c cc ccc cccc 3 2 1 3 2 1 333 22322 1131211 000 00 0 cnnxn = dn ® xn = dn cnn cn-1,n-1xn-1 +cn-1,nxn = dn-1 ® xn-1 = dn-1 - cn-1,nxn cn-1,n-1 Introdução – Sistemas Triangulares Sistema Triangular Superior , ù é = ù é ´ ù é nnnn n n n d d d d x x x x c cc ccc cccc 3 2 1 3 2 1 333 22322 1131211 000 00 0 c22x2 + c23x3 + ...+c2nxn = d2 ® x2 = d2 - c23x1 -...- c2nxn c22 c11x1 + c12x2 +c13x3 + ...+c1nxn = d1 ® x1 = d1 -c12 -c13x1 -...-c1nxn c11 Introdução – Sistemas Triangulares Sistema Triangular Superior As substituições retroativas podem ser representadas por: .1...,,1,, 1 -= - = += nni c xcd x ii n ij jiji i ù é = ù é ´ ù é nnnn n n n d d d d x x x x c cc ccc cccc 3 2 1 3 2 1 333 22322 1131211 000 00 0 Introdução – Sistemas Triangulares Sistema Triangular Superior Exemplo: Determinar a solução do sistema triangular superior utilizando as substituições retroativas: , , ù é - = ù é ù é - - 8 28 2 1 2000 5400 4730 1625 4 3 2 1 x x x x 2x4 = 8® x4 = 8 2 = 4 4x3 + 5x4 = 28® x3 = 28- 5´ 4 4 = 2 3x2 + 7x3 - 4x4 = -2 ® x2 = -2- 7´ 2+ 4´ 4 3 = 0 18/02/2019 6 Introdução – Sistemas Triangulares Sistema Triangular Superior Exemplo: Determinar a solução do sistema triangular superior utilizando as substituições retroativas: Logo, o vetor solução é dado por: ù é - = ù é ù é - - 8 28 2 1 2000 5400 4730 1625 4 3 2 1 x x x x 3 5 426021 1625 14321 -= -´-´+ =®=++- xxxxx x = -3 0 2 4 é ù Introdução – Métodos de Solução Os métodos numéricos para a solução de sistemas lineares podem ser divididos em dois grupos: Métodos Diretos Fornecem a solução do sistema, caso ela exista, após um número finito de iterações (soluções arredondadas também podem ocorrer) Métodos Iterativos Geram uma sequência de vetores {x(k)} a partir de uma aproximação inicial x(0). Sob certas condições esta sequência converge para a solução x do sistema, caso ela exista Métodos Diretos Pertencem a essa classe todos métodos utilizados no primeiro e segundo graus Esses métodos não são eficientes para a resolução de sistemas lineares de grande porte, ou seja, sistemas que envolvam um grande número de equações e variáveis Para o caso de sistemas lineares n x n, com solução única, o vetor x é dado por x = A-1b, onde A-1 é a inversa da matriz de coeficientes A. O cálculo de A-1 é demorado e, por isso, não competitivo com os métodos que veremos a seguir: Eliminação de Gauss e Decomposição LU Eliminação de Gauss Consiste em transformar o sistema linear original em um sistema linear triangular superior equivalente Resolução do novo sistema utilizando as substituições retroativas A solução encontrada para o sistema equivalente será a mesma do sistema linear original Conceito de Sistemas Equivalentes Eliminação de Gauss A transformação do sistema linear original em outro equivalente é feita através das seguintes operações elementares: Trocar duas equações Multiplicar uma equação por uma constante não nula Adicionar um múltiplo de uma equação a uma outra equação î í ì -=- =+ ® î í ì =+ -=- 222 94 94 222 21 21 21 21 xx xx xx xx î í ì =+ -=- ® î í ì =+ -=- 94 1 94 222 21 21 21 21 xx xx xx xx î í ì -=- =+ ® î í ì -=- =+ 1 832 1 94 21 21 21 21 xx xx xx xx Eliminação de Gauss - Execução Passo 1: Construção da matriz aumentada Ab Importância: É necessário transformar matriz A em uma matriz triangular superior Todavia, todas as operações elementares aplicadas sobre as linhas de A, também devem ser refletidas no vetor de termos independentes b ù é = nnnnnn n n baaaa baaa baaa Ab 321 222221 111211 18/02/2019 7 Eliminação de Gauss - Execução Passo 2: Eliminar os coeficientes de x1 presentes nas linhas 2,3,...,n, fazendo assim a21 = a31, = ... = an1 = 0, sendo a11 chamado de pivô e a linha 1 de linha pivotal Substituir a linha 2, L2, pela combinação linear Substituir a linha 3, L3, pela combinação linear: 11 21 2112122 :, a a mqualnaLmLL =×-= 11 31 3113133 :, a a mqualnaLmLL =×-= Eliminação de Gauss - Execução Passo 2: Continuar a substituição até a linha n Caso algum elemento app=0, achar outra linha k onde akp≠ 0 e trocar tais linhas. Caso a linha k não exista, o sistema linear não possui solução Próximos Passos: Eliminar os coeficientes de x2 nas linhas 3, 4, ..., n (fazendo a32=a42=...=an2 = 0) Eliminar os coeficientes de x3 nas linhas 4, 5, ..., n (fazendo a43=a53=...=an3 = 0) e assim sucessivamente Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: Matriz aumentada: ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx ù é -- - - = 1132 3344 5132 Ab Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: Pivô da linha 1: 2 ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx 7120513223344 2 2 11 21 2112122 ---=-´--= ==´-= L a a mLmLL 6260513211132 1, 3 11 31 3113133 --=-´---= ==´-= L a a mLmLL Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: Obtemos então a seguinte matriz aumentada: ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx ù é -- --- - = 6260 7120 5132 Ab Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: Pivô da linha 2: -2 ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx 15500712036260 3, 3 22 32 3223233 =---´---= ==´-= L a a mLmLL 18/02/2019 8 Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: Nova matriz [Ab] e sistema linear equivalente obtido: ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx ï î ï í ì = -=-- =-+ 155x 7x2x 5x3x2x 3 32 321 ù é --- - = 15500 7120 5132 Ab Eliminação de Gauss Exemplo: Resolver o sistema linear abaixo: O novo sistema obtido é resolvido utilizando-se as substituições retroativas: Logo, o vetor solução é dado por: ï î ï í ì -=+- =-+ =-+ 132 3344 532 321 321 321 xxx xxx xxx 3x155x 33 =®= ù é = 3 2 1 x 1x22x5362x5x3x2x 111321 =®=®=-+®=-+ 2x42x732x7x2x 22232 =®-=-®-=--®-=-- Eliminação de Gauss No método de Gauss os multiplicadores das linhas são gerados a partir da seguinte fórmula: sendo aii o pivô e aik o elemento a ser zerado Assim, podemos concluir: O método de Gauss não funciona quando o pivô é nulo Quando o pivô é muito próximo de zero, os multiplicadores gerados para as linhas são muito grandes, ocasionando um aumento nos erros de arredondamento gerados durante a execução do método. Solução: Pivoteamento Parcial nikni a a m ii ik ik ...,,1...,,1 +=== Eliminação de Gauss - Pivoteamento Parcial Melhoria do Método de Gauss Consiste em escolher o elemento de maior valor (em módulo) em cada coluna para ser o pivô Garante que os multiplicadores estarão sempre entre 0 e 1 Minimiza a amplificação de erros de arredondamento durante as eliminações Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: Matriz aumentada: ï î ï í ì =+- -=-+- =+- 29564 15182 1123 321 321 321 xxx xxx xxx ù é - --- - = 29564 15182 11231 Ab Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: Maior elemento (em módulo) da primeira coluna: 4. Logo este será o primeiro pivô. Assim: 75,375,05,102956425,011231 25,0 1 31 11 1331311 -=-´--= ==´-= L a a mLmLL ù é - --- - = 29564 15182 11231 Ab L2 = L2 -m23 ´ L3, m23 = a21 a31 = -0, 5 L2 = -2 8 -1 -1 5 é ù - (-0,5)´ 4 -6 5 29 é ù = 0 5 1,5 -0, 5 é ù 18/02/2019 9 Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: Obtemos então a seguinte matriz aumentada: ù é - - - = 29564 5,05,150 75,375,05,10 Ab ï î ï í ì =+- -=-+- =+- 29564 15182 1123 321 321 321 xxx xxx xxx Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: Maior elemento (em módulo) da segunda coluna: 5. Logo este será o segundo pivô. Assim: L1 = L1 -m12 ´ L2 m12 = a12 a22 = -0,3 L1 = 0 -1,5 0, 75 3, 75 é ù - (-0,3)´ 0 5 1, 5 -0, 5 é ù = 0 0 1, 2 3,6 é ù ù é - - - = 29564 5,05,150 75,375,05,10 Ab Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: Nova matriz [Ab]: Trocando a ordem das linhas, chegamos ao seguinte sistema equivalente: ï î ï í ì =+- -=-+- =+- 29564 15182 1123 321 321 321 xxx xxx xxx ù é - -= 29564 5,05,150 6,32,100 Ab ï î ï í ì = -=+ =+- 6,32,1 5,05,15 29564 3 32 321 x xx xxx Eliminação de Gauss - Pivoteamento Parcial Exemplo: Resolver o sistema linear abaixo: O novo sistema obtido é resolvido utilizando-se as substituições retroativas: Logo, o vetor solução é dado por: ù é -= 3 1 2 x x1 -3x2 + 2x3 =11 -2x1 +8x2 -1x3 = -15 4x1 - 6x2 + 5x3 = 29 ì í ï î ï 1, 2x3 = 3, 6® x3 = 3 5x2 +1, 5x3 = -0, 5® 5x2 + 4, 5 = -0, 5® 5x2 = -5® x2 = -1 4x1 - 6x2 + 5x3 = 29® 4x1 + 6+15= 29 ® 4x1 = 8® x1 = 2 Exercício Resolva o sistema linear abaixo utilizando o método direto de Eliminação de Gauss com pivoteamento parcial. ï î ï í ì =++ =++ =++ 5734 22 8425 321 321 321 xxx xxx xxx Determinante O determinante da matriz de coeficientes pode ser obtido através da matriz triangular resultante da aplicação da Eliminação de Gauss. Basta considerar no cálculo a influência das operações elementares realizadas durante o processo de eliminação Vamos então analisar essas relações: 1) Se duas linhas de uma matriz A forem trocadas, então o determinante da nova matriz B será: )det()det( AB -= 10)det( 22 41 10)det( 41 22 -=® ù é - ==® ù é - = BBeAA 18/02/2019 10 Determinante 2) Se todos os elementos de uma linha de A forem multiplicados por uma constante k, então o determinante da matriz resultante B será: 3) Se um múltiplo escalar de uma linha de A for somado a outra linha, então o determinante da nova matriz B será: det(B) = k´det(A) A= 1 4 2 -2 é ù ® det(A) = -10 e B= 1 4 1 -1 é ù ® det(B) = -5 det(B) = det(A) Determinante 4) Se A for uma matriz triangular ou diagonal de ordem n, então seu determinante será igual ao produto dos elementos da diagonal principal, ou seja: = == n i iinn aaaaaA 1 332211 ...)det( A= 2 3 0 -1 é ù ® det(A) = -2 e B= 3 0 0 0 -5 0 0 0 -1 é ù ® det(B) =15 Determinante 5) Se uma matriz A for multiplicada por uma matriz B, o determinante da matriz resultante C será: )det()det()det( BAC ´= 3)det( 11 03 10)det( 43 21 =® ù é ==® ù é - = BBeAA 30)det( 413 21 =® ù é - = CC Determinante Exemplo: Calcular o determinante da matriz utilizada no último exemplo: Matriz de coeficientes: Depois de 3 combinações lineares das linhas e uma troca de linhas, chegamos à seguinte matriz triangular: ù é - -- - = 564 182 231 A ï î ï í ì -=+- -=-+- =+- 29564 15182 1123 321 321 321 xxx xxx xxx ù é - = 2.100 5.150 564 A Determinante Exemplo: Calcular o determinante da matriz utilizada no último exemplo: Pela propriedade 3, não há alteração no determinante de B, todavia, pela propriedade 1, , assim:)det()det( AB -= 24)2,154()det()det( -=´´-=-= BA Exercício Resolva o sistema linear abaixo utilizando o método direto de Eliminação de Gauss. Use a técnica de pivoteamento parcial se necessário. Solução: ï ï î ï ï í ì =+++ =-++ =+++ =+++ 7239335 181284 25154193 8426 4321 4321 4321 4321 xxxx xxxx xxxx xxxx ù é- = 1 11 20 138 x 18/02/2019 11 Sistemas Lineares Decomposição LU Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Decomposição LU O objetivo é fatorar a matriz dos coeficientes A em um produto de duas matrizes L e U. Seja: A = matriz de coeficientes do sistema linear ù é = nnnnn n n aaaa aaa aaa A 321 22221 11211 Decomposição LU e o produto LU: sendo: L = matriz triangular inferior unitária U = matriz triangular superior ù é ´ ù é = nn n n n nnn u uu uuu uuuu lll ll l LU 000 00 0 1 0 01 001 0001 333 22322 1131211 321 3231 21 ( )ilii = ,1 Decomposição LU tem-se então: Logo, o sistema Ax = b pode ser reescrito como Ax = b LUx = b Fazendo Ux = y, a equação acima reduz-se a Ly = b. Resolvendo o sistema triangular inferior (utilizando as substituições sucessivas) Ly = b, obtém-se o vetor y ù é ´ ù é == ù é = nn n n n nnn nnnnn n n u uu uuu uuuu lll ll l LU aaaa aaa aaa A 000 00 0 1 0 01 001 0001 333 22322 1131211 321 3231 21 321 22221 11211 Decomposição LU O vetor y é então utilizado como termo independente no sistema triangular superior Ux = y, cuja solução x é calculada pelas substituições retroativas A Decomposição LU é um dos processos mais empregados. Uma das vantagens é que podemos resolver qualquer sistema linear que tenha A como matriz de coeficientes. Se o vetor b for alterado, a solução do novo sistema linear será quase que imediata Decomposição LU - Execução Exemplo: Resolver o sistema abaixo, utilizando a Decomposição LU: Passo 1: Aplicar o método da Eliminação de Gauss à matriz de A. Pivô linha 1: 1 ù é -= ù é ´ ù é - -- - 29 15 11 564 182 231 3 2 1 x x x 18/02/2019 12 Decomposição LU - Execução Exemplo: Obtemos então a seguinte matriz de coeficientes: Pivô linha 2: 2 3602314564 4, 3 11 31 3113133 -=-´--= ==´-= L a a mLmLL ù é - - 360 320 231 12003203360 3 3 22 32 3223233 -=´--= ==´-= L a a mLmLL Decomposição LU - Execução Exemplo: Nova matriz de coeficientes: A matriz L é então constituída pelos multiplicadores utilizados nas eliminações de cada uma das linhas, logo: ù é - - 1200 320 231 ù é -® ù é ® ù é = 134 012 001 1 01 001 1 01 001 3231 21 3231 21 mm m ll lL Decomposição LU - Execução Exemplo: Nova matriz de coeficientes: A matriz U é própria matriz de coeficientes, obtida após a Eliminação de Gauss: ù é - - 1200 320 231 ù é - - = ù é = 1200 320 231 00 0 33 2322 131211 u uu uuu U Decomposição LU - Execução Exemplo: Assim: Substituindo a matriz de coeficientes A no sistema, temos então LUx = b. Fazendo Ux = y, temos então Ly = b. Assim o próximo passo na solução do problema é calcular o valor do vetor y. ù é - - ´ ù é -= ù é - -- - ®= 1200 320 231 134 012 001 564 182 231 LUA Decomposição LU - Execução Exemplo: Passo 2: Calcular a solução do sistema Ly = b: Logo: ù é -= ù é ´ ù é - 29 15 11 134 012 001 3 2 1 y y y 111 =y 711215152 2221 =®´+-=®-=+- yyyy 3673114292934 33321 -=®´-´-=®=++ yyyyy Ty 36711 -= Decomposição LU - Execução Exemplo: Passo 3: De posse do valor de y, calcular então a solução do sistema Ux = y: Logo: ù é - = ù é ´ ù é - - 36 7 11 1200 320 231 3 2 1 x x x 33612 33 =®-=- xx 12/)337(732 2232 -=®´-=®=+ xxxx 232)1(3111123 11321 =®´--´+=®=+- xxxxx Tx 312 -= 18/02/2019 13 Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Solução: ï î ï í ì =++ =++ =++ 3234 22 1423 321 321 321 xxx xxx xxx ù é- = 0 5 3 x Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Demonstrando a solução: Aplicando a Eliminação de Gauss à matriz A: Pivô linha 1: 3 ï î ï í ì =++ =++ =++ 3234 22 1423 321 321 321 xxx xxx xxx ù é = 234 211 423 A 3/23/104233/1211 3 1 2 11 21 2112122 =´-= ==´-= L a a mLmLL Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Obtemos então a seguinte matriz de coeficientes: Pivô linha 2: 1/3 3/103/104233/4234 3 4 , 3 11 31 3113133 -=´-= ==´-= L a a mLmLL ù é - = 3/103/10 3/23/10 423 A Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Nova matriz de coeficientes: A matriz L é então constituída pelos multiplicadores utilizados nas eliminações de cada uma das linhas, logo: ù é - = 400 3/23/10 423 A ù é ® ù é ® ù é = 113/4 013/1 001 1 01 001 1 01 001 3231 21 3231 21 mm m ll lL Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Nova matriz de coeficientes: A matriz U é própria matriz de coeficientes, obtida após a Eliminação de Gauss: ù é - = 400 3/23/10 423 A ù é - = ù é = 400 3/23/10 423 00 0 33 2322 131211 u uu uuu U Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Assim: Substituindo a matriz de coeficientes A no sistema, temos então LUx = b. Fazendo Ux = y, temos então Ly = b. Assim o próximo passo na solução do problema é calcular o valor do vetor y. ù é - ´ ù é = ù é ®= 400 3/23/10 423 113/4 013/1 001 234 211 423 LUA 18/02/2019 14 Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Calcular a solução do sistema Ly = b: Fazendo os cálculos, vamos encontrar: ù é = ù é ´ ù é 3 2 1 113/4 013/1 001 3 2 1 y y y Ty 03/51= Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. De posse do valor de y, calcular então a solução do sistema Ux = y: Fazendo os cálculos, vamos encontrar: Tx 053-= ù é = ù é ´ ù é - 0 3/5 1 400 3/23/10 423 3 2 1 x x x Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Mas e se, ao invés de colocar o termo de multiplicação na linha pivotal, eu colocá-lo na própria linha cujo elemento desejo zerar? Funciona? ï î ï í ì =++ =++ =++ 3234 22 1423 321 321 321 xxx xxx xxx ù é = 234 211 423 A Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Mas e se, ao invés de colocar o termo de multiplicação na linha pivotal, eu colocá-lo na própria linha cujo elemento desejo zerar? Funciona? Aplicando a Eliminação de Gauss à matriz A: Pivô linha 1: 3 2104232113 3 2 2112212 --=+´-= =+-= L mLLmL 4/104/10423234 4 3 4 3 , 3 3113313 -=+-= =´-= L mLLmL Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Mas e se, ao invés de colocar o termo de multiplicação na linha pivotal, eu colocá-lo na própria linha cujo elemento desejo zerar? Funciona? Obtemos então a seguinte matriz de coeficientes: Pivô linha 2: -1 ù é - --= 4/104/10 210 423 A 12002104/104/104 4 3 3223323 -=--+-´-= =+-= L mLLmL Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Mas e se, ao invés de colocar o termo de multiplicação na linha pivotal, eu colocá-lo na própria linha cujo elemento desejo zerar? Funciona? Nova matriz de coeficientes: ù é - --= 1200 210 423 A 18/02/2019 15 Exercício Resolva o sistema linear abaixo utilizando o método direto da Decomposição LU. Mas e se, ao invés de colocar o termo de multiplicação na linha pivotal, eu colocá-lo na própria linha cujo elemento desejo zerar? Funciona? Fazendo então A = LU Percebemos, pela igualdade acima, que, ao fazer L x U, não encontramos A. Daí concluímos que a multiplicação deve sempre ser feita na linha pivotal, como mostra a fórmula trabalhada em sala. ù é - --´ ù é -- -= ù é ®= 1200 210 423 144/3 013 001 234 211 423 LUA Decomposição LU – Pivoteamento Parcial Os motivos para Pivoteamento Parcial na Decomposição LU são os mesmos de sua utilização na Eliminação de Gauss: Evitar pivô nulo Evitar que os multiplicadores mij tenham valores muito grandes Decomposição LU – Pivoteamento Parcial No pivoteamento parcial, a decomposição é feita da forma: PA = LU onde P é uma matriz de permutações que será construída das linhas de uma matriz identidade I, colocadas na mesma ordem das linhas que geram a matriz triangular superior U. A matriz L é formada pelos multiplicadores utilizados na eliminação nas respectivas linhas de U. Assim, para resolver o sistema Ax = b, tem-se: Ax = b PAx = Pb LUx = Pb Fazendo Ux = y, então Ly = Pb Decomposição LU – Pivoteamento Parcial Exemplo: Resolver o sistema abaixo, utilizando a Decomposição LU, com Pivoteamento Parcial: Passo 1: Aplicar o método da Eliminação de Gauss, com Pivoteamento Parcial, à matriz A. Primeiro pivô: 4 ù é -= ù é ´ ù é - -- - 29 15 11 564 182 231 3 2 1 x x x Decomposição LU - Pivoteamento Parcial Exemplo: Obtemos então a seguinte matriz de coeficientes: Segundo pivô: 5 5,150564)5,0(182 5,0, 2 31 11 2332322 =-´----= -==´-= L a a mLmLL ù é - - 564 5,150 75,05,10 2,1005,150)3,0(75,05,10 3,0 1 32 12 1221211 =´---= -==´-= L a a mLmLL Decomposição LU - Pivoteamento Parcial Exemplo: Nova matriz de coeficientes: A matriz L é então constituída pelos multiplicadores relativos a cada uma das linhas pivotais, logo: ù é - ® ù é - 2,100 5,150 564 564 5,150 2,100 ù é - -® ù é ® ù é = 13,025,0 015,0 001 1 01 001 1 01 001 1213 23 3231 21 mm m ll lL 18/02/2019 16 Decomposição LU - Pivoteamento Parcial Exemplo: Nova matriz de coeficientes: A matriz U é própria matriz de coeficientes, obtida após o pivoteamento: ù é - ® ù é - 2,100 5,150 564 564 5,150 2,100 ù é - = ù é = 2,100 5,150 564 00 0 33 2322 131211 u uu uuu U Decomposição LU - Pivoteamento Parcial Exemplo: Nova matriz de coeficientes: A matriz P possui as linhas de uma matriz identidade na ordem das linhas pivotais. P pode ser vista ainda como uma matriz similar à identidade com as linhas colocadas de modo que os elementos iguais a 1 estejam nas colunas relativas aos índices das linhas pivotais. ù é - ® ù é - 2,100 5,150 564 564 5,150 2,100 Decomposição LU - Pivoteamento Parcial Exemplo: Assim: Assim, para resolver o sistema Ax = b, temos: Ax = b PAx = Pb LUx = Pb Fazendo Ux = y, então Ly = Pb ù é - ù é - -= ù é - -- - ù é ®= 2,100 5,150 564 13,025,0 015,0 001 564 182 231 001 010 100 LUPA Decomposição LU - Pivoteamento Parcial Exemplo: Passo 2: Calcular a solução do sistema Ly = Pb: A multiplicação Pb ordena as linhas de b na ordem das linhas pivotais Logo: ù é -= ù é ´ ù é - - 11 15 29 13,025,0 015,0 001 3 2 1 y y y 291 =y 5,0295,015155,0 2221 -=®´+-=®-=+- yyyy 5,03,02925,011113,025,0 3321 -´+´-=®=+- yyyy Ty 6,35,029 -= 6,33 =y Decomposição LU - Pivoteamento Parcial Exemplo: Passo 3: De posse do valor de y, calcular então a solução do sistema Ux = y: Logo: ù é -= ù é ´ ù é - 6,3 5,0 29 2,100 5,150 564 3 2 1 x x x 36,32,1 33 =®= xx 15/)35,15,0(5,05,15 2232 -=®´--=®-=+ xxxx ( ) 4/35)1(62929564 1321 ´--´+=®=+- xxxx Tx 312 -= 21 =x Decomposição LU - Pivoteamento Parcial Exemplo: Passo 4 (Opcional): Verificação da exatidão do resultado obtido através do vetor resíduo r = b – Ax: Logo, a solução x obtida é exata. A verificação acima pode ser utilizada também para validar os resultados encontrados pelos outros métodos estudados até agora. ù é = ù é -´ ù é - -- - - ù é -®-= 0 0 0 3 1 2 564 182 231 29 15 11 Axbr 18/02/2019 17 Determinante Considerando que: PA = LU det (PA) = det (LU) pela propriedade dos determinantes vista anteriormente: matriz triangular produto dos pivôs troca de linhas necessárias para transformar a matriz de permutações P em uma matriz identidade. )det( )det()det( )det( P UL A = 1...)det( 1 332211 === = n i iinn lllllL = = n i iiuU 1 )det( tP )1()det( -= Determinante Considerando que: PA = LU det (PA) = det (LU) pela propriedade dos determinantes vista anteriormente: Logo: )det( )det()det( )det( P UL A = = = -= - ´ = n i ii t t n i ii u u A 1 1 )1( )1( 1 )det( Determinante Exemplo: Calcular o determinante da matriz utilizada no último exemplo: Para calcular o determinante, precisamos encontrar o valor de t, isto é, o número de trocas de linhas necessárias para transformar a matriz P em uma matriz identidade. Voltando na matriz, percebemos que somente uma troca é suficiente. Assim, temos t=1 e: 242,154)1()1()det( 1 1 -=´´´-=-= = n i ii t uA ù é - ® ù é - -- - = 2,100 5,150 564 564 182 231 A Sistemas com Matriz Singular Quando a matriz de coeficientes do sistema linear for singular, ou seja, det(A) = 0, o sistema pode ter infinitas soluções ou não ter solução. Será mostrado como diferenciar essas situações. Exemplo: Resolver os sistemas Ax = b e Ax = c utilizando a decomposição LU com pivoteamento parcial, sendo ù é -= ù é -= ù é - -- - = 80 10 20 10 12 22 , 151 182 231 cebA Sistemas com Matriz Singular Exemplo: Os três fatores são: Para Ax=b, a solução do sistema Ly = Pb é dada por: ù é = ù é -- = ù é -= 100 001 010 000 5,110 182 , 115,0 015,0 001 PeUL ù é- =® ù é- = ù é ´ ù é - 0 16 12 10 22 12 115,0 015,0 001 3 2 1 y y y y Sistemas com Matriz Singular Exemplo: A solução do sistema Ux = y é dada por: Logo: ù é- = ù é ´ ù é -- 0 16 12 000 5,110 182 3 2 1 x x x )(00 333 soluçãoéxdevalorqualquerxx =®= 5,116165,1 232 -=®=+ xxx ( ) 2/)5,116(8121282 1321 -+---=®-=-+- xxxx 5,6701 -=x 18/02/2019 18 Sistemas com Matriz Singular Exemplo: Assim, o vetor solução do sistema é dado por , ou seja, o sistema Ax=b apresenta infinitas soluções, uma para cada valor de . Para resolver o sistema Ax=c, não é necessário calcular novamente L, U e P. Como a matriz de coeficientes A é comum aos dois sistemas (Ax=b e Ax=c), os cálculos feitos anteriormente podem ser reaproveitados. Tx 5,1165,670 --= | Sistemas com Matriz Singular Exemplo: Assim, para Ax = c, solução de Ly = Pc é ù é- =® ù é- = ù é ´ ù é - 70 15 10 80 20 10 115,0 015,0 001 3 2 1 y y y y Sistemas com Matriz Singular Exemplo: A solução do sistema Ux = y é dada por: Logo: Assim, o sistema Ax = c não tem solução pois tal que . ù é- = ù é ´ ù é -- 70 15 10 000 5,110 182 3 2 1 x x x xxx ®®= 33 700 3x 00 3 x Exercício Resolva o sistema linear pela Decomposição LU, utilizando o pivoteamento parcial, e verificar a exatidão e unicidade da solução: Solução: ï î ï í ì -=- =++ =+- 234 322 943 31 321 321 xx xxx xxx ù é -= 2 1 1 x Sistemas Lineares Métodos Iterativos Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Métodos Iterativos A solução de problemas complexos com sistemas lineares tende à geração/existência de matrizes de coeficientes grandes e/ou esparsas Grandes Comum para n > 100.000 Esparsas Maioria dos coeficientes nulos Resolução de sistemas esparsos por métodos diretos Processos de triangularização e fatoração Onerosos, por não preservarem a esparsidade original, que pode ser útil por facilitar a resolução do sistema. 18/02/2019 19 Métodos Iterativos Métodos mais apropriados para a resolução de sistemas de natureza esparsa Métodos Iterativos Gauss-Jacobi, Gauss-Seidel Métodos Iterativos Os métodos iterativos consistem em gerar, a partir de um vetor inicial x0, uma sequência de vetores {x0, x1, x2, …, xk, …} que deve convergir para a solução x do sistema ù é )0( )0( 3 )0( 2 )0( 1 nx x x x ù é )1( )1( 3 )1( 2 )1( 1 nx x x x ù é )2( )2( 3 )2( 2 )2( 1 nx x x x ù é )( )( 3 )( 2 )( 1 k n k k k x x x x Métodos Iterativos Lembretes importantes: Como todo processo iterativo, estes métodos sempre apresentarão um resultado aproximado, que será tão próximo do resultado real conforme o número de iterações realizadas. Além disto, também é preciso ter cuidado com a convergência destes métodos. Métodos Iterativos - Funcionamento Os métodos iterativos funcionam a partir da transformação do sistema linear Ax = b em x = Cx + g, onde: A: matriz dos coeficientes (n x n) x: vetor das variáveis (n x 1) b: vetor dos termos constantes, (n x 1) C: matriz n x n g: vetor n x 1 Métodos Iterativos - Funcionamento Conhecida a estimativa inicial, x(0), obtém-se consecutivamente os vetores: De um modo geral, a aproximação x(k+1) é calculada pela fórmula: x(k+1) = C x(k) + g, k = 0, 1, ... chamada de função de iteração, dada na forma matricial o)aproximaçã ésima-(k , o)aproximaçã (segunda , o)aproximaçã (primeira , )1()( )1()2( )0()1( gCxx gCxx gCxx kk += += += - Sistemas Lineares Gauss - Jacobi 18/02/2019 20 Método de Gauss - Jacobi Dado o sistema linear: e supondo , i = 1, …, n. ï ï î ï ï í ì =+++ =+++ =+++ nnnnnn nn nn bxaxaxa bxaxaxa bxaxaxa ... ... ... 2211 22222121 11212111 0iia Método de Gauss - Jacobi Isolamos então o vetor x mediante a separação pela diagonal. Assim, a partir da primeira equação do sistema: obtemos: e, analogamente: 11212111 ... bxaxaxa nn =+++ )...( 1 13132121 11 1 nn xaxaxab a x +---= )...( 1 23231212 22 2 nn xaxaxab a x ---= )...( 1 113311 -----= nnnnnn nn n xaxaxab a x Método de Gauss - Jacobi Dessa forma, temos xk+1 = C xk + g, onde: x(k+1) C x(k) g ù é + ù é ù é --- --- --- --- = ù é nnnnnnnnnnnnn n n n n ab ab ab ab x x x x aaaaaa aaaaaa aaaaaa aaaaaa x x x x / / / / 0 0 0 0 333 222 111 3 2 1 321 33333323331 22222232221 11111131112 3 2 1 Método de Gauss - Jacobi O método de Gauss-Jacobi consiste em, dado , aproximação inicial, obter através da relação recursiva : 0x ( ) ( )......,0 kxx ( ) ( ) gCxx kk +=+1 ï ï ï ï î ï ï ï ï í ì +---= +---= +---= -- + + + )...( 1 )...( 1 )...( 1 )( 11, )( 22 )( 11 )1( )( 2 )( 323 )( 1212 22 )1( 2 )( 1 )( 313 )( 2121 11 )1( 1 k nnn k n k nn nn k n k nn kkk k nn kkk xaxaxab a x xaxaxab a x xaxaxab a x Método de Gauss - Jacobi O processo é repetido até que o vetor esteja suficientemente próximo ao vetor A distância entre duas iterações é dada por assim, dada uma precisão , o vetor será escolhido como , solução aproximada da solução exata, se Podemos utilizar também como critério de parada o erro relativo: ( )kx ( )1-kx - max d 1)(k-i (k) i (k) xx= ( )kx x d(k) e e max d d )( (k) (k) r e = k ix Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi com : x0 = 0,7 -1,6 0,6 é ù e e = 0, 05 ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx 18/02/2019 21 Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : O processo iterativo é dado por: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ( ) ( )ï ï ï î ï ï ï í ì ++--=--= --+-=---= +--=--= + + + 10 6 0 10 3 10 2 3 26 10 1 5 8 5 1 0 5 1 8 5 1 10 7 10 1 10 2 0) 2 (7 10 1 )( 3 )( 2 )( 1 )( 2 )( 1 )1( 3 )( 3 )( 2 )( 1 )( 3 )( 1 )1( 2 )( 3 )( 2 )( 1 )( 3 )( 2 )1( 1 kkkkkk kkkkkk kkkkkk xxxxxx xxxxxx xxxxxx Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Na forma matricial temos: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ù é = 03/10– 1/5- 1/5-01/5- 1/10 -2/10 -0 C ( ) ( ) gCxx kk +=+1 ù é = 6/10 8/5- 7/10 ge Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Assim para k=0 e temos: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ï î ï í ì =+-´-´-=+--= -=-´-´-=---= =+´--´-=+--= 94,06,0)6,1(3,07,02,0 6,0 3,02,0 86,16,16,02,07,02,06,1 2,0 2,0 96,07,06,0 1,0 )6,1(2,07,0 1,0 2,0 )0( 2 )0( 1 )1( 3 )0( 3 )0( 1 )1( 2 )0( 3 )0( 2 )1( 1 xxx xxx xxx x 0 = 0,7 -1,6 0,6 é ù Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Logo: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ù é =+= 0,94 1,86- 0,96 g C (0)(1) xx e=== 1828,0 86,1 34,0 max 34,0 d )1( )1( i r x Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Calculando , para temos: d(1)r ù é = ù é = 0,94 1,86- 0,96 6,0 6,1- 7,0 (1))0( xex ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx 34,0 - 26,0 - 26,0 - )0( 3 )1( 3 )0( 2 )1( 2 )0( 1 )1( 1 = = = xx xx xx Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Assim para k=1 e temos: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ï î ï í ì -=+-´-´-=+--= -=-´-´-=---= =+´--´-=+--= 966,06,0)86,1(3,096,02,0 6,0 3,02,0 98,16,194,02,096,02,06,1 2,0 2,0 978,07,094,0 1,0 )86,1(2,07,0 1,0 2,0 )1( 2 )1( 1 )2( 3 )1( 3 )1( 1 )2( 2 )1( 3 )1( 2 )2( 1 xxx xxx xxx x1 = 0, 96 -1,86 0, 94 é ù 18/02/2019 22 Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Logo: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ù é =+= 0,966 1,98- 0,978 g C )1()2( xx e=== 0,0606 98,1 12,0 max 12,0 d )2( )2( i r x Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Calculando , para : d(2)r ù é = ù é = 0,966 1,98- 0,978 0,94 1,86- 0,96 )2()1( xex ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx 026,0 - 12,0 - 018,0 - )1( 3 )2( 3 )1( 2 )2( 2 )1( 1 )2( 1 = = = xx xx xx Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Prosseguindo com as iterações temos: Para k=2 ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx e 0,0163 1,9888 0,0324 d 0,9984 1,9888- 0,9994 g C (2) r (2)(3) == ù é =+= xx Método de Gauss - Jacobi Exemplo: Resolva o sistema abaixo utilizando Gauss- Jacobi : Logo a solução obtida pelo método de Gauss-Jacobi é: ï î ï í ì =++ -=++ =++ 6 10 3 2 8 5 7 2 10 321 321 321 xxx xxx xxx ù é == 0,9984 1,9888- 0,9994 (3)xx Método de Gauss – Jacobi - Convergência No exemplo estudado, o valor de foi fornecido como entrada do problema. Todavia, a convergência ou não dos métodos iterativos independe da aproximação inicial escolhida Teorema: Critério das Linhas Dado um sistema Ax=b, é condição suficiente para a convergência do método iterativo de Gauss-Jacobi: ou seja, o somatório do módulo de todos os elementos da linha, exceto o elemento da diagonal principal, deve ser menor que este elemento (0)x n ..., 3, 2, 1,i para , 1 = = ii n ij j ij aa Analisando a matriz A do sistema linear do exemplo anterior: Assim: Logo, como temos a convergência garantida para o método de Gauss-Jacobi Método de Gauss – Jacobi - Convergência ù é = 1032 151 1210 A 53210 2115 31210 323133 232122 131211 =+=+= =+=+= =+=+= aaa aaa aaa 3 2, 1,i para 1 = = ii n ij j ij aa 18/02/2019 23 Exemplo: Dado o sistema: O critério das linhas não é satisfeito pois: Contudo, se permutarmos a primeira equação com a segunda, temos o sistema linear: Método de Gauss – Jacobi - Convergência ï î ï í ì -=+ =++ -=++ 686 3225 23 32 321 321 xx xxx xxx 4131 131211 =+=+= aaa ï î ï í ì -=+ -=++ =++ 686 23 3225 32 321 321 xx xxx xxx Exemplo: Dado o sistema: O novo sistema é equivalente ao sistema original e sua matriz A satisfaz o critério de linhas: Método de Gauss – Jacobi - Convergência ï î ï í ì -=+ =++ -=++ 686 3225 23 32 321 321 xx xxx xxx ù é = 860 131 225 A Conclusão: Sempre que o critério de linhas não for satisfeito, devemos tentar uma permutação de linhas e/ou colunas de forma a obtermos uma disposição para a qual a matriz dos coeficientes satisfaça o critério de linhas Método de Gauss – Jacobi - Convergência Calcule as 3 primeiras iterações do método de Gauss- Jacobi do sistema linear abaixo: Utilize como chute inicial: Exercício ï î ï í ì =+- =-- =-+ 5152 3865 7924 321 321 321 xxx xxx xxx ù é = 3 2 1 )0(x Calcule as 3 primeiras iterações do método de Gauss- Jacobi do sistema linear abaixo: Utilize como chute inicial: Resp: Exercício ï î ï í ì =+- =-- =-+ 5152 3865 7924 321 321 321 xxx xxx xxx ù é = 3 2 1 )0(x ù é -= 533,0 667,3 5,7 )1(x Sistemas Lineares Gauss - Seidel 18/02/2019 24 Método de Gauss-Seidel Similarmente ao método de Gauss-Jacobi, conhecida a estimativa inicial, x(0), obtém-se consecutivamente os vetores x(1), x(2), ..., x(k) Todavia, ao se calcular xj (k+1), usa-se todos os valores x1 (k+1), x2 (k+1), ..., xj-1 (k+1) que já foram calculados e os valores xj+1 (k), xj+2 (k), ..., xn (k) restantes. Método de Gauss-Seidel O processo do método de Gauss Seidel se dá a partir das equações: ( ) ( ) ( ) ( ) ï ï ï ï ï ï î ï ï ï ï ï ï í ì ----= -----= -----= -----= + -- +++ -- +++ -- ++ -- + 1 11, 1 22 1 11 1 311,3 1 232 1 1313 33 1 3 211,2323 1 1212 22 1 2 111,13132121 11 1 1 ... 1 ... 1 ... 1 ... 1 k nnn k n k nn nn k n k nn k nn kkk k nn k nn kkk k nn k nn kkk xaxaxab a x xaxaxaxab a x xaxaxaxab a x xaxaxaxab a x Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel com : 05,0 0 0 0 0 = ù é = eex ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: O processo iterativo é dado por: ( ) ( )ï ï ï î ï ï ï í ì --=--= --=--= --=--= +++++ +++ + 6 3 6 3 6 0 3 30 6 1 4 1 4 3 4 6 36 4 1 5 1 5 1 5 5 ) 5( 5 1 )1( 2 )1( 1 )1( 2 )1( 1 )1( 3 )( 3 )1( 1 )( 3 )1( 1 )1( 2 )( 3 )( 2 )( 3 )( 2 )1( 1 kkkkk kkkkk kkkkk xxxxx xxxxx xxxxx ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Assim para k=0 e temos: ï î ï í ì -=´-´-=--= =´-´-=--= =´-´-=--= 875,075,05,015,00 5,05,00 75,0025,0175,05,1 25,0 75,05,1 10 2,0 02,01 2,0 2,01 )1( 2 )1( 1 )1( 3 )0( 3 )1( 1 )1( 2 )0( 3 )0( 2 )1( 1 xxx xxx xxx ù é = 0 0 0 0x ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Logo: ù é - =+= 875,0 75,0 1 g C (0)(1) xx ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx 18/02/2019 25 875,0 - 75,0 - 1 - )0( 3 )1( 3 )0( 2 )1( 2 )0( 1 )1( 1 = = = xx xx xx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Calculando , para temos: d(1)r ù é - = ù é = 875,0 75,0 1 0 0 0 (1))0( xex ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx e=== 1 1 1 max 1 d )1( )1( i r x Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Assim para k=1 e temos: ï î ï í ì -=´-´-=--= =-´-´-=--= =-´-´-=--= 9875,095,05,0025,15,00 5,05,00 95,0875,025,0025,175,05,1 25,0 75,05,1 025,1875,0 2,0 75,02,01 2,0 2,01 )2( 2 )2( 1 )2( 3 )1( 3 )2( 1 )2( 2 )1( 3 )1( 2 )2( 1 xxx xxx xxx ù é = 0,875- 0,75 1 )1(x ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Logo: ù é - =+= 9875,0 95,0 025,1 g C )1()2( xx ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx e=== 1951,0 025,1 2,0 max 2,0 d )2( )1( i r x Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Calculando , para : d(2)r ù é - = ù é - = 9875,0 95,0 025,1 875,0 75,0 1 )2()1( xex ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx 1125,0 - 20,0 - 025,0 - )1( 3 )2( 3 )1( 2 )2( 2 )1( 1 )2( 1 = = = xx xx xx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Prosseguindo com as iterações temos: Para k=3 ε 4090,0 d 9930,9 9912,0 0075,1 g C )3( r )2()3( = ù é - =+= xx ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx Método de Gauss - Seidel Exemplo: Resolva o sistema abaixo utilizando Gauss- Seidel: Prosseguindo com as iterações temos: Logo a solução obtida pelo método de Gauss-Seidel é: ï î ï í ì =++ =++ =++ 0 6 3 3 6 4 3 5 5 321 321 321 xxx xxx xxx ù é - == 9930,9 9912,0 0075,1 )3(xx 18/02/2019 26 Método de Gauss – Seidel - Convergência Para o método de Gauss Seidel, utilizaremos os seguintes critérios de convergência Critério de Sassenfeld Critérios das Linhas Critério de Sassenfeld Sejam o valores dados por: n - ordem do sistema linear que se deseja resolver aij - coeficientes das equações do sistema Este critério garante que o método de Gauss-Seidel convergirá para um dado SEL se o valor M, definido por: for menor que 1 (M<1). Além disso, quanto menor o valor de mais rápida será a convergência n ..., 3, 2, i para 11 1 1 12 1 11 1 = ù é +´´=´= += - == n ij ij i j jij ii i n j j aa a ea a bbb ib i M ni bmax 1 = b Critério de Sassenfeld Seja A a matriz dos coeficientes e b o vetor dos termos constantes, dados por: ù é 444434241 334333231 224232221 114131211 baaaa baaaa baaaa baaaa b1 = 1 a11 × | a12 |+ | a13 |+ | a14 |( ) b2 = 1 a22 × a21 b1 + a23 + a24( ) b3 = 1 a33 × a31 b1 + a32 b2 + a34( ) b4 = 1 a44 × a41 b1 + a42 b2 + a43 b3( ) Critério de Sassenfeld Exemplo: Mostrar que a solução do sistema a seguir convergirá pelo método de Gauss-Seidel. ï ï î ï ï í ì -=+++- =++-- -=--+ =+-+ 10 48,0 2,1 4,0 1 2,0 2,0 1,0 8,7 3,06,0 3 6,0 4,0 2,02,0 2 4321 4321 4321 4321 xxxx xxxx xxxx xxxx Critério de Sassenfeld Exemplo: Mostrar que a solução do sistema a seguir convergirá pelo método de Gauss-Seidel. ( ) ( ) ( ) ( ) 2736,0358,08,044,02,17,04,0 4 1 358,02,044,02,07,01,0 1 1 44,03,06,07,06,0 3 1 7,02,02,01 2 1 4 3 2 1 =´+´+´´= =+´+´´= =++´´= =++×= b b b b 4,0 0,8 1,2 0,4 0,2 1,0 0,2- 0,1- 0,3- 0,6- 3,0 0,6 0,2 0,2- 1,0 2,0 A = 7,0max 1 == ini M b Logo, como M < 1, o método de Gauss – Seidel converge para o sistema em questão Critério de Sassenfeld O exemplo anterior também satisfaz o critério das linhas: 4,0 0,8 1,2 0,4 0,2 1,0 0,2- 0,1- 0,3- 0,6- 3,0 0,6 0,2 0,2- 1,0 2,0 A = 4,28,02,14,04 5,02,02,01,01 5,13,06,06,03 4,12,02,012 43424144 34323133 24232122 14131211 =++=++= =++=++= =++=++= =++=++= aaaa aaaa aaaa aaaa 4 3, 2, 1,i para 1 = = ii n ij j ij aa 18/02/2019 27 Considerações finais – Gauss Seidel Tanto o Critério de Sassenfeld quanto o Critério das Linhas são condições suficientes, porém não necessárias, para a convergência do método de Gauss-Seidel para um dado sistema linear Um dado sistema pode não satisfazer estes critérios e ainda convergir Um sistema pode não satisfazer o Critério das Linhas, porém sua convergência será garantida se satisfizer o Critério de Sassenfeld Considerações finais – Gauss Seidel Critério das Linhas x Critério de Sassenfeld Exemplo: Seja o sistema linear abaixo: O Critério das Linhas não é satisfeito, visto que: Todavia o Critério de Sassenfeld é satisfeito, uma vez que: Convergência garantida! î í ì =+ =+ 18 2 6 23 10 21 21 xx xx 62 2122 == aa ( ) 3,01,06 2 1 1,01 10 1 21 =´´==´= bb e Considerações finais – Métodos Diretos/Iterativos Métodos Diretos Processos finitos Teoricamente obtêm a solução de qualquer sistema linear ( det(A) diferente de 0 ) Podem sofrer com problemas de arredondamento Solução: Técnicas de Pivoteamento Métodos Iterativos Convergem para a solução do sistema linear somente sob certas condições Sofrem menos com problemas de arredondamento Convergência independe do valor de x(0) Somente erros cometidos na última iteração afetam a solução Calcule as 3 primeiras iterações do método de Gauss- Siedel do sistema linear abaixo: Utilize como chute inicial: Resp. Exercício ï î ï í ì =+- =-- =-+ 5152 3865 7924 321 321 321 xxx xxx xxx ù é = 3 2 1 )0(x ù é = 0,067 75,1 5,7 )1(x ù é = 3,0 266,0 026,1 )2(x ù é = 315,0 01,1 292,2 )3(x 05_Interpolacao_Folheto.pdf 18/02/2019 1 Interpolação Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Programa 1. Introdução 2. Métodos de Interpolação a) Linear b) Quadrática c) Lagrange d) Newton e) Gregory-Newton 3. Estudo do Erro Introdução Muitas funções são conhecidas apenas em um conjunto finito e discreto de pontos de um intervalo [a,b] Exemplo: A tabela abaixo informa o número de carros que passam por um determinado pedágio em um determinado dia A partir dos dados acima, suponhamos que se queira calcular o número de carros que passariam pelo pedágio às 11:30 Horário 10:00 11:00 12:00 13:00 14:00 15:00 Quantidade Carros (x1000) 2,69 1,64 1,09 1,04 1,49 2,44 Introdução A interpolação tem o objetivo de ajudar na resolução deste tipo de problema Interpolar uma função y = f(x), em um conjunto discreto de pontos, pertencentes a um intervalo [a, b], consiste em substituí-la, ou aproximá-la, por uma outra função, y = g(x), escolhida dentro de uma classe de funções definida a priori e que satisfaça algumas propriedades. A função y = g(x) é, então, usada no lugar da função y = f(x) Introdução A necessidade de se efetuar esta substituição surge em várias situações, como por exemplo: Quando y = f(x) não é conhecida na sua forma analítica, nas apenas em um conjunto discreto de pontos (xi, yi), i = 0, 1, ..., n Quando y = f(x) é conhecida na sua forma analítica, mas operações como a diferenciação e a integração são difíceis (ou impossíveis) de realizar, ou seja, a função é difícil de ser tratada. Logo, podemos procurar uma outra função que seja uma aproximação y = f(x) e cujo manuseio seja bem mais simples Introdução Teoricamente, a função interpoladora, y = g(x), pode ser de tipos variados, tais como: polinomiais trigonométricas exponenciais logarítmicas Porém, para consideraremos apenas o estudo das funções polinomiais 18/02/2019 2 Conceito Formal de Interpolação Numérica Consideremos (n+1) pontos distintos: x0, x1, ..., xn, chamados nós da interpolação, e os valores de f(x) nesses pontos: f(x0), f(x1), ..., f(xn) Uma forma de interpolação de f(x) consiste em se obter uma determinada função g(x) tal que: )()( )()( )()( )()( 22 11 00 nn xfxg xfxg xfxg xfxg = = = = Conceito Formal de Interpolação Numérica Graficamente, para n=5 temos (6 pontos), temos: Nos nós da interpolação as funções f(x) e g(x) assumem os mesmos valores! Interpolação Polinomial Dados os pontos (x0, f(x0)), (x1, f(x1)), ..., (xn, f(xn)), portanto (n+1) pontos, queremos aproximar f(x) por um polinômio pn(x) de grau menor ou igual a n, tal que: f(xk) = pn(xk) k = 0,1, 2 , …. , n Representaremos pn(x) por: Portanto, obter pn(x) significa obter os coeficientes a0, a1, ..., an n nn xaxaxaaxp ++++= ...)( 2 210 Interpolação Polinomial Surgem então as seguintes perguntas: Existe sempre um polinômio que satisfaça as condições anteriores? Caso exista, ele é único? Teorema: Existe um único polinômio pn(x) de grau ≤ n, tal que pn(xk) = f(xk), k = 0,1, 2 , …. , n desde que xk ≠ xj, j ≠ k Demonstrando o Teorema: Seja Da condição , montamos o seguinte sistema linear: com n+1 equações e n+1 incógnitas a0, a1, ..., an Interpolação Polinomial nkxfxp kkn ...,,2,1,0),()( == a0 +a1x0 +a2x0 2 +...+ anx0 n = f (x0 ) a0 +a1x1 +a2x1 2 +...+ anx1 n = f (x1) a0 +a1xn +a2xn 2 +...+ anxn n = f (xn) ì í n nn xaxaxaaxp ++++= ...)( 2 210 Interpolação Polinomial Demonstrando o Teorema: Da notação Matricial temos: v x a = f onde V é uma matriz de Vardermonde e, portanto, desde que x0, x1, x2, …, xn sejam pontos distintos, temos det(v) ≠ 0. Logo, o sistema linear admite solução única. A matriz coluna a é a matriz das incógnitas e a matriz coluna f é a matriz das constantes f(xi) = yi = = = )( )( )( , 1 1 1 1 0 1 0 2 1 2 11 0 2 00 nn n nnn n n xf xf xf fe a a a a xxx xxx xxx v 18/02/2019 3 Interpolação Polinomial Existem várias formas de se obter o polinômio pn(x) Linear Quadrática Lagrange Newton Gregory-Newton, etc Interpolação Interpolação Linear Interpolação Linear Dados dois pontos distintos de uma função y = f(x) : (x0, y0) e (x1, y1), deseja-se calcular o valor de ӯ para um determinado valor de x entre x0 e x1, utilizando-se a interpolação polinomial O polinômio interpolador é uma unidade menor que o número de pontos conhecidos Assim, nesse caso, o polinômio interpolador terá grau 1, ou seja: p1(x) = a0 + a1x Interpolação Linear Para determinar este polinômio, os coeficientes a0 e a1 devem ser calculados de forma que se tenha: p1(x0) = f(x0) = y0 e p1(x1) = f(x1) = y1 Ou seja, basta resolver o sistema linear abaixo: onde a1 e a0 são as incógnitas Reescrevendo o sistema na forma matricial, temos: í ì =+ =+ 1110 0010 yxaa yxaa = 1 0 1 0 1 0 1 1 y y a a x x Transformando em um sistema triangular equivalente (utilizando a Eliminação de Gauss), temos: A solução do sistema é dada por: 0101 00 0 1 yyxx yx Interpolação Linear 11 00 1 1 yx yx í ì = =+ 01011 0010 yyxxa yxaa 0100 01 01 1 xayae xx yy a = = Interpolação Linear Logo, o polinômio interpolador pode ser escrito da seguinte forma: ou, de forma mais apropriada: xaaxp 101 )( += )()( 0 01 01 01 xx xx yy yxp += )( 010 xxay +=xaxay 1010 )( += 18/02/2019 4 Interpolação Linear O determinante da matriz A é diferente de zero, sempre que x0 ≠ x1, logo para pontos distintos o sistema tem solução única O polinômio interpolador p1(x) = a1x+ a0 tem como imagem geométrica uma reta, portanto a função f(x) está sendo aproximada por uma reta que passa pelos pontos conhecidos (x0, y0) e (x1, y1) Interpolação Linear O gráfico abaixo, mostra geometricamente, os dois pontos (x0, y0) e (x1, y1), e a reta que passa por eles. Interpolação Linear Exemplo: Seja a função y = p1(x) definida pelos pontos da tabela abaixo, determinar o valor de p1(15): Solução: Aplicando a fórmula temos: Logo: i 0 1 x 10 20 y = f(x) 250 432 )()( 0 01 01 01 xx xx yy yxp += )10( 1020 250432 250)(1 += xxp 34168152,18)15(1 =+=p 682,18 += x Interpolação Interpolação Quadrática Interpolação Quadrática Se conhecermos três pontos distintos de uma função, então o polinômio interpolador será: p2(x) = a2x 2 + a1x + a0 O polinômio p2(x) é conhecido como função quadrática cuja imagem geométrica é uma parábola Portanto, a função f(x) é aproximada por uma parábola que passa pelos três pontos conhecidos (x0, y0), (x1, y1) e (x2, y2) Interpolação Quadrática Para determinar os valores de a0, a1 e a2 é necessário resolver o sistema: onde a1, a0 e a2 são as incógnitas e os pontos (x0, y0), (x1, y1) e (x2, y2) são conhecidos. í ì =++ =++ =++ 2 2 22210 1 2 12110 0 2 02010 yxaxaa yxaxaa yxaxaa 18/02/2019 5 Interpolação Quadrática Reescrevendo o sistema na forma matricial, temos: O sistema acima admite solução única, pois det(X) = (x2 –x0) (x2 – x1)(x1 – x0) ≠ 0. Logo, pelos três pontos ordenados passa um único polinômio interpolador de 2° grau = 2 1 0 2 1 0 2 22 2 11 2 00 1 1 1 y y y a a a xx xx xx Interpolação Quadrática Para encontrar os coeficientes ai, basta resolver o sistema de equações ( Eliminação de Gauss, Decomposição LU, Método de Gauss-Jacobi, Método de Gauss-Siedel, etc ) Interpolação Quadrática Exemplo: Seja a função y = p2(x) definida pelos pontos da tabela abaixo, determinar o valor de p2(0,2): Solução: Montando o sistema linear: i 0 1 2 x 0,1 0,6 0,8 y = f(x) 1,221 3,320 4,953 í ì =++ =++ =++ 953,48,08,0 320,36,06,0 221,11,01,0 2 210 2 210 2 210 aaa aaa aaa Interpolação Quadrática Exemplo: Seja a função y = p2(x) definida pelos pontos da tabela abaixo, determinar o valor de p2(0,2): Solução: Montando o sistema linear: = 953,4 320,3 221,1 8,08,01 6,06,01 1,01,01 2 1 0 2 2 2 a a a i 0 1 2 x 0,1 0,6 0,8 y = f(x) 1,221 3,320 4,953 Interpolação Quadrática Exemplo: Seja a função y = p2(x) definida pelos pontos da tabela abaixo, determinar o valor de p2(0,2): Solução: Aplicando a Eliminação de Gauss ao sistema, chegaremos ao seguinte sistema triangular superior: í ì = =+ =++ 567,01,0 733,363,07,0 221,101,01,0 567,01,000 733,363,07,00 221,101,01,01 2 21 210 a aa aaa i 0 1 2 x 0,1 0,6 0,8 y = f(x) 1,221 3,320 4,953 Interpolação Quadrática Exemplo: Seja a função y = p2(x) definida pelos pontos da tabela abaixo, determinar o valor de p2(0,2): Solução: Resolvendo o sistema triangular superior, encontramos: a2 = 5,667; a1 = 0,231; a0 = 1,141 Logo o polinômio terá a seguinte forma: 2 2 667,5231,0141,1)( xxxp ++= i 0 1 2 x 0,1 0,6 0,8 y = f(x) 1,221 3,320 4,953 18/02/2019 6 Interpolação Quadrática Exemplo: Seja a função y = p2(x) definida pelos pontos da tabela abaixo, determinar o valor de p2(0,2): Solução: Calculando o valor de p2(0,2) 2 2 2,0667,52,0231,0141,1)2,0( ++=p 414,1)2,0(2 =p i 0 1 2 x 0,1 0,6 0,8 y = f(x) 1,221 3,320 4,953 Exercício prático Seja a tabela: Usando interpolação quadrática sobre pontos adequados, calcule f(0.35), onde f(x) = ex. x 0 0.1 0.2 0.3 0.4 0.5 ex 1 1.11 1.22 1.35 1.49 1.65 Interpolação Interpolação de Lagrange Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Interpolação de Lagrange As interpolações vistas anteriormente são casos particulares da interpolação de Lagrange. Vamos estudar agora o polinômio interpolador de grau menor ou igual a n, sendo dados n + 1 pontos distintos Interpolação de Lagrange Sejam x0, x1, x2, ..., xn, ( n + 1 ) pontos distintos e yi = f(xi), i = 0, 1, ..., n Seja pn(x) o polinômio de grau ≤ n que interpola f em x0, ..., xn. Podemos representar pn(x) na forma: onde os polinômios Lk(x) são de grau n. )(...)()()( 1100 xLyxLyxLyxp nnn +++= Interpolação de Lagrange Para cada i, queremos que a condição pn(xi) = yi seja satisfeita, ou seja: A forma mais simples de se satisfazer esta condição é impor: í ì = = ikse ikse xL ik 1 0 )( iinniiin yxLyxLyxLyxp =+++= )(...)()()( 1100 18/02/2019 7 Interpolação de Lagrange E para isso, definimos Lk(x) por: como o numerador de Lk(x) é um produto de n fatores da forma: (x - xi), i = 0, 1, 2, ..., n, i ≠ k, então Lk(x) é um polinômio de grau n e, assim, pn(x) é um polinômio de grau menor ou igual a n Note que, na fórmula acima, e como pré-definido anteriormente ))...()()...()(( ))...()()...()(( )( 1110 1110 nkkkkkkk nkk k xxxxxxxxxx xxxxxxxxxx xL = + + 1)( =kk xL kise0)( =ik xL Interpolação de Lagrange Assim, para x = xi, i = 0, …, n temos: Então a interpolação de Lagrange para o polinômio interpolador é: Fazendo a substituição de Lk(x) na fórmula, temos: = === n k iiiiikkin yxLyxLyxp 0 )()()( = = == n k n kjj jk j kkkn xx xx xLondexLyxp 0 ,0 )( )( )(),()( = = = n k n kjj jk j kn xx xx yxp 0 ,0 )( )( )( Interpolação de Lagrange Exemplo: Considere a tabela de dados abaixo: Solução: Pela forma de Lagrange, temos que: onde: i 0 1 2 x -1 0 2 y = f(x) 4 1 -1 )()()()( 2211002 xLyxLyxLyxp ++= ))(( ))(( )( 2010 21 0 xxxx xxxx xL = 3 22 xx = )21)(01( )2)(0( = xx Interpolação de Lagrange Exemplo: Considere a tabela de dados abaixo: Solução: Pela forma de Lagrange, temos que: i 0 1 2 x -1 0 2 y = f(x) 4 1 -1 ))(( ))(( )( 2101 20 1 xxxx xxxx xL = ))(( ))(( )( 1202 10 2 xxxx xxxx xL = 2 22 = xx )20)(10( )2)(1( + + = xx 6 2 xx + = )02)(12( )0)(1( + + = xx Interpolação de Lagrange Exemplo: Considere a tabela de dados abaixo: Solução: Voltando na forma de Lagrange, temos: Chegamos então a: i 0 1 2 x -1 0 2 y = f(x) 4 1 -1 + + + = 6 )1( 2 2 1 3 2 4)( 222 2 xxxxxx xp )()()()( 2211002 xLyxLyxLyxp ++= 2 2 3 2 3 7 1)( xxxp += Interpolação de Lagrange Exercício: Considere a tabela de dados abaixo: Escreva o polinômio de Lagrange de ordem 2 para esse conjunto de pontos. Calcule p2(2,5). Resposta: p2(x) = 20,897 – 23,497x + 12,354x 2 p2(2,5) = 39,37 i 0 1 2 x 1,1 2,2 3,5 y = f(x) 10 29 90 18/02/2019 8 Interpolação Linear/Quadrática Exercício: Dados os pontos abaixo, encontre: a) O polinômio interpolador de ordem 2 (Parábola) que ajusta os pontos acima, utilizando o método de Gauss para triangularizar o sistema de equações Resposta: p2(x) = -0,5775 + 0,7567x – 0,0214x 2 b) O polinômio de ordem 1 ( reta ) que ajusta os pontos x0 e x2 Resposta: p1(x) =0,7059 + 0,2647x i 0 1 2 x 3 9 20 y = f(x) 1,5 4,5 6,0 Interpolação Interpolação de Newton Elaborado pelo Prof. Wellington Passos de Paula Adaptado por Marconi de Arruda Pereira Interpolação de Newton Um polinômio pn(x) com grau n é baseado em n + 1 pontos Um novo ponto é observado. Como atualizar pn(x) passando a pn+1(x)? Não queremos recalcular todos os valores... É possível fazer isto? Sim: Interpolação de Newton Interpolação de Newton A forma de Newton para o polinônio pn(x), que interpola f(x) em (n+1) pontos distintos x0, x1, ..., xn , é a seguinte: No Método da Interpolação de Newton, os valores de dk são dados por Diferenças Divididas de ordem k. 110 102010 ....... .......)( + ++++= nn n xxxxxxd xxxxdxxddxp Seja f(x) definida em (n+1) pontos distintos x0, x1, ..., xn. O operador Diferenças Divididas é dado: Dizemos que f[x0, x1, x2, ..., xn] é a diferença dividida de ordem k da função f(x) sobre os k + 1 pontos Interpolação de Newton – Dif. Divididas )(][ 00 xfxf = 01 01 01 01 10 )()(][][ ],[ xx xfxf xx xfxf xxf = = 02 1021 210 ],[],[ ],,[ xx xxfxxf xxxf = 0 121021 210 ],....,,,[],...,,[ ],...,,,[ xx xxxxfxxxf xxxxf n nn n = Interpolação de Newton – Dif. Divididas Conhecidos os valores que f(x) assume nos pontos distintos x0, x1, x2, ..., xn, podemos construir a seguinte tabela das Diferenças Divididas: )(][ 00 xfxf = x Ordem 0 Ordem 1 Ordem 2 Ordem n x0 x1 x2 ..... ...... ...... ......... xn ],[ 10 xxf ],[ 21 xxf ],,[ 210 xxxf ][ 0xf ][ 1xf ][ 2xf ],[ 32 xxf ],[ 1 nn xxf ],,[ 321 xxxf ],..,,[ 10 nxxxf 18/02/2019 9 Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 1 : calcule f[x0, x1] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] x2 f [x2] x3 f [x3] Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 2 : calcule f[x1, x2] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] f [x1, x2] x2 f [x2] x3 f [x3] Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 3 : calcule f[x2, x3] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] f [x1, x2] x2 f [x2] f [x2, x3] x3 f [x3] Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 4 : calcule f[x0, x1, x2] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] f [x0, x1 ,x2] f [x1, x2] x2 f [x2] f [x2, x3] x3 f [x3] Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 5 : calcule f[x1, x2, x3] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] f [x0, x1 ,x2] f [x1, x2] x2 f [x2] f [x1, x2 ,x3] f [x2, x3] x3 f [x3] Interpolação de Newton – Dif. Divididas Montando a tabela das Diferenças Divididas para n = 4: Passo 6 : calcule f[x0, x1, x2, x3] x Ordem 0 Ordem 1 Ordem 2 Ordem 3 x0 f [x0] f [x0, x1] x1 f [x1] f [x0, x1 ,x2] f [x1, x2] f [x0, x1 ,x2, x3] x2 f [x2] f [x1, x2 ,x3] f [x2, x3] x3 f [x3] 18/02/2019 10 Interpolação de Newton – Dif. Divididas Propriedades das Diferenças Divididas: Satisfazem a propriedade de ser simétrica aos argumentos: Exemplo: ],[ ][][][][ ],[ 01 10 10 01 01 10 xxf xx xfxf xx xfxf xxf = = = ].......,,[],,[],,[ 021201210 xxxfxxxfxxxf == Interpolação de Newton – Dif. Divididas Propriedades das Diferenças Divididas: Cada coeficiente dn do Polinômio Interpolador de Newton corresponde ao operador de grau n de Diferenças Divididas: nn dxxxxf dxxxf dxxf dxf = = = = ],...,,,[ ],,[ ],[ ][ 210 2210 110 00 Interpolação de Newton Logo a forma de Newton para o polinômio de ordem n que interpola f(x) é dada por: ))....()(](,..,,,[... ...))(](,,[)](,[][)( 110210 102100100 ++ +++= nn n xxxxxxxxxxf xxxxxxxfxxxxfxfxp 110 102010 ....... .......)( + ++++= nn n xxxxxxd xxxxdxxddxp Interpolação de Newton Exemplo: Usando a forma de Newton, calcule o polinômio p4(x) que interpola f(x) nos pontos abaixo: Tabela das Diferenças Divididas: x -1 0 1 2 3 f(x) 1 1 0 -1 -2 Interpolação de Newton Exemplo: Usando a forma de Newton, calcule o polinômio p4(x) que interpola f(x) nos pontos abaixo: A forma de Newton que interpola estes pontos é dada por: x -1 0 1 2 3 f(x) 1 1 0 -1 -2 )2)(1)(0)(1( 24 1 )1)(0)(1( 6 1 )0)(1( 2 1 )1(01)(4 + + + + + +++= xxxx xxx xxxxp Interpolação de Newton Exemplo: Usando a forma de Newton, calcule o polinômio p4(x) que interpola f(x) nos pontos abaixo: A forma de Newton que interpola estes pontos é dada por: x -1 0 1 2 3 f(x) 1 1 0 -1 -2 )2)(1)(1( 24 1 )1)(1( 6 1 )1( 2 1 1)(4 + + + + + += xxxx xxx xxxp 18/02/2019 11 Interpolação de Newton Lagrange e Newton acham o mesmo polinômio, apenas expressam e calculam de forma diferente O cálculo de Newton e estável numericamente (como Lagrange) mas é computacionalmente mais eficiente que Lagrange (menos operações) Interpolação de Newton Exercício: Usando a forma de Newton, calcule o polinômio p2(x) que interpola f(x) nos pontos abaixo: Resposta: i 0 1 2 x -1 0 2 y = f(x) 4 1 -1 2 2 3 2 3 7 1)( xxxp += Interpolação Interpolação de Gregory-Newton Interpolação de Gregory-Newton Muitas vezes são encontrados problemas de interpolação cuja tabela de pontos conhecidos tem valores que são igualmente espaçados, ou seja: Assim xi+1 – xi = h, para todo i, sendo h uma constante hxxxxxxxx nn ===== 1231201 ... hxx ii += 1 hixxi *0 += Interpolação de Gregory-Newton Diferenças Ordinárias ou Finitas )()()(1 xfhxfxf += )()(0 xfxf = )()()( 112 xfhxfxf += )()()( 11 xfhxfxf nnn += Interpolação de Gregory-Newton Diferenças Ordinárias ou Finitas x Ordem 0 Ordem 1 Ordem 2 … Ordem n x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) x2 f (x2) ∆ 2 f (x1) ... ... ... ... ... ∆n f (x0) xn-2 f (xn-2) ∆1 f (xn-2) xn-1 f (xn-1) ∆ 2 f (xn-2) ∆1 f (xn-1) xn f (xn) 18/02/2019 12 Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 1: calcule ∆1 f (x0) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) x2 f (x2) x3 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 2: calcule ∆1 f (x1) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆1 f (x1) x2 f (x2) x3 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 3: calcule ∆1 f (x2) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆1 f (x1) x2 f (x2) ∆1 f (x2) x3 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 4: calcule ∆1 f (x3) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆1 f (x1) x2 f (x2) ∆1 f (x2) x3 f (x3) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 5: calcule ∆2 f (x0) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) x2 f (x2) ∆1 f (x2) x3 f (x3) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 6: calcule ∆2 f (x1) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) x2 f (x2) ∆ 2 f (x1) ∆1 f (x2) x3 f (x3) ∆1 f (x3) x4 f (x4) 18/02/2019 13 Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 7: calcule ∆2 f (x2) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) x2 f (x2) ∆ 2 f (x1) ∆1 f (x2) x3 f (x3) ∆ 2 f (x2) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 8: calcule ∆3 f (x0) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) ∆ 3 f (x0) x2 f (x2) ∆ 2 f (x1) ∆1 f (x2) x3 f (x3) ∆ 2 f (x2) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 9: calcule ∆3 f (x1) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) ∆ 3 f (x0) x2 f (x2) ∆ 2 f (x1) ∆1 f (x2) ∆ 3 f (x1) x3 f (x3) ∆ 2 f (x2) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Montando a tabela das Diferenças Ordinárias ou Finitas para n = 4: Passo 10: calcule ∆4 f (x0) x Ordem 0 Ordem 1 Ordem 2 Ordem 3 Ordem 4 x0 f (x0) ∆1 f (x0) x1 f (x1) ∆ 2 f (x0) ∆1 f (x1) ∆ 3 f (x0) x2 f (x2) ∆ 2 f (x1) ∆ 4 f (x0) ∆1 f (x2) ∆ 3 f (x1) x3 f (x3) ∆ 2 f (x2) ∆1 f (x3) x4 f (x4) Interpolação de Gregory-Newton Relação entre Diferenças Divididas e Diferenças Ordinárias Teorema: Se xj = x0 + j * h, para j = 0, 1, 2, …, n, então: n n n hn xf xxxxf ! )( ...,,,, 0210 = Interpolação de Gregory-Newton Relação entre Diferenças Divididas e Diferenças Ordinárias Prova: por indução, podemos mostrar que esta regra é válida para valores maiores que 2 )( 00 xfxf = = = 02 1021 210 ],[],[ ,, xx xxfxxf xxxf = = = 01 01 01 01 10 )()(][][ , xx xfxf xx xfxf xxf h xf h xfhxf )()()( 000 = + = 2 0 2 01 2 )( 2 )()( h xf h h xf h xf = = 18/02/2019 14 Interpolação de Gregory-Newton Partindo da fórmula original do Método de Newton, que é: podemos derivar a fórmula do Método de Gregory- Newton, que utiliza as Diferenças Ordinárias: ))....()(](,..,,,[ ...))(](,,[)](,[][)( 110210 102100100 + +++= nn n xxxxxxxxxxf xxxxxxxfxxxxfxfxp ))....()(( ! )( ...))(( 2 )( )( )( )()( 110 0 102 0 2 0 0 0 + + + += nn n n xxxxxx hn xf xxxx h xf xx h xf xfxp Interpolação de Gregory-Newton Exemplo: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p4(0,5) pela fórmula de Gregory-Newton: Solução: x -1 0 1 2 3 f(x) 2 1 2 5 10 Exemplo: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p4(0,5) pela fórmula de Gregory-Newton: Solução: A forma de Gregory-Newton que interpola estes pontos é dada por: Interpolação de Gregory-Newton x -1 0 1 2 3 f(x) 2 1 2 5 10 ))()()(( !4 )( ))()(( !3 )( ))(( !2 )( )( !1 )( )()( 32104 0 4 2103 0 3 102 0 2 01 0 04 xxxxxxxx h xf xxxxxx h xf xxxx h xf xx h xf xfxp + + + + += Exemplo: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p4(0,5) pela fórmula de Gregory-Newton: Solução: Interpolação de Gregory-Newton x -1 0 1 2 3 f(x) 2 1 2 5 10 )2)(1)()(1( 124 0 )1)()(1( 16 0 ))(1( 12 2 )1( 11 1 2)(4 + ++ + + ++ += xxxxxxx xxxxp 1)( 24 += xxp Exemplo: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p4(0,5) pela fórmula de Gregory-Newton: Solução: Interpolação de Gregory-Newton x -1 0 1 2 3 f(x) 2 1 2 5 10 25,1125,01)5,0( 24 =+=+= xp Grau do Polinômio Interpolador Para a escolha do grau do polinômio interpolador: 1) Construa a tabela de diferenças divididas 2) Examine as diferenças na vizinhança do ponto de interesse Se as diferenças de ordem k forem praticamente constantes, ou se as diferenças de ordem k+1 variarem em torno de zero, o polinômio de grau k será o que melhor aproximará a função na região considerada 18/02/2019 15 Grau do Polinômio Interpolador Exemplo: Dada com os valores tabelados: Um polinômio de grau 1 é uma boa aproximação para Vamos provar? xxf =)( x 1 1,01 1,02 1,03 1,04 1,05 y = f(x) 1 1,005 1,01 1,0149 1,0198 1,0247 xxf =)( Grau do Polinômio Interpolador a) Tabela de diferenças: Exercício prático Seja a tabela: Usando interpolação de Gregory Newton, gere os polinômios P4(x) e P2(x) sobre pontos adequados, calcule f(0.35), P4(0.35) e P2(0.35) onde f(x) = e x, e calcule o erro obtido pelas funções interpoladoras. x 0 0.1 0.2 0.3 0.4 0.5 ex 1 1.11 1.22 1.35 1.49 1.65 Interpolação Estudo do Erro Estudo do Erro na Interpolação Ao se aproximar uma função f(x) por um polinômio interpolador pn(x), de grau n, comete-se um erro En (x) tal que seu valor estimado é: onde é a derivada de ordem (n + 1) da função interpolada e )!1( )( ))...()(()()()( )1( 10 + == + n f xxxxxxxpxfxE x n nnn ),( 0 nx xx )()1( x nf + Estudo do Erro na Interpolação Interpolação linear de f1(x) e f2(x) x f(x) x0 x1 f1(x0)= f2(x0)=p1(x0) f1(x) p1(x) f2(x)f1(x1)= f2(x1)=p1(x1) 18/02/2019 16 Estudo do Erro na Interpolação Interpolação linear de f1(x) e f2(x) O mesmo polinômio p1(x) interpola f1(x) e f2(x) em x0 e x1 O erro E1 1(x)=f1(x) - p1(x) > E1 2(x)= f2(x) - p1(x) para todo x de (x0 , x1) O erro depende da concavidade da curva, ou seja, de f1”(x) e f2”(x) Estudo do Erro na Interpolação Exemplo: Considere a função e a tabela dada abaixo. Obtenha p1(0,7) por interpolação linear e calcule o erro cometido: Solução: Pela fórmula de Newton temos: 1)( += xexf x x 0 0,5 1 1,5 2,0 y = f(x) 0 1,1487 2,7183 4,9811 8,3890 ],[)()()( 10001 xxfxxxfxp += Estudo do Erro na Interpolação Exemplo: Considere a função e a tabela dada abaixo. Obtenha p1(0,7) por interpolação linear e calcule o erro cometido: Solução: Como x = 0,7, logo x0 = 0,5 e x1 = 1, temos: 1)( += xexf x x 0 0,5 1 1,5 2,0 y = f(x) 0 1,1487 2,7183 4,9811 8,3890 += 5,01 1487,17183,2 )5,0(1487,1)(1 xxp 1392,3)5,0(1487,1)(1 += xxp 7765,1=1392,3)5,07,0(1487,1)7,0(1 +=p Estudo do Erro na Interpolação Exemplo: Considere a função e a tabela dada abaixo. Obtenha p1(0,7) por interpolação linear e calcule o erro cometido: Solução: O cálculo do erro então será: 1)( += xexf x x 0 0,5 1 1,5 2,0 y = f(x) 0 1,1487 2,7183 4,9811 8,3890 == )7,0()7,0()7,0( 11 pfE == 7765,17137,1 0628,0= 0628,00628,0)7,0(1 ==E Estudo do Erro na Interpolação A fórmula do erro tem uso limitado na prática pois raras são as situações em que conhecemos . Além disso, o ponto nunca é conhecido Assim o que normalmente se usa é uma fórmula para a estimativa do erro, que é dada pela expressão abaixo: onde )()1( xf n+ x !1 )).....()(()()()( 110 + = + n M xxxxxxxpxfxE nnnn .1max )!1( 1 += + + nordemdivididasdiferenças n Mn Estudo do Erro na Interpolação Exemplo: Sejam os valores de f(x) informados na tabela abaixo: a) Utilize a fórmula de Newton para obter f(0,47) usando um polinômio de grau 2 b) Encontre uma estimativa para o erro x 0,2 0,34 0,4 0,52 0,6 0,72 y = f(x) 0,16 0,22 0,27 0,29 0,32 0,37 18/02/2019 17 Estudo do Erro na Interpolação a) Gerando a tabela de Diferenças Divididas: Estudo do Erro na Interpolação a) Escolhendo Substituindo: 6,0,52,0,4,0 210 === xxx ],,[))((],[)()()( 2101010002 xxxfxxxxxxfxxxfxp ++= )04115,1()52,0)(4,0()1667,0()4,0(27,0)(2 ++= xxxxp )47,0(2780,0)47,0(2 fp = Estudo do Erro na Interpolação b) Tabela de Diferenças Divididas: Estudo do Erro na Interpolação b) como Logo |2492,18||)6,047,0)(52,047,0)(4,047,0(||)47,0(| 2 E 3 2 10303,8|)47.0(| E 009,0278,0)47,0(2 =p !12 ))()(()( 122102 + + M xxxxxxxE .1max )!1( 1 += + + nordemdivididasdiferenças n Mn Exercício prático Seja a tabela: Usando interpolação de Gregory Newton, gere os polinômios P4(x) e P2(x) sobre pontos adequados, calcule f(0.35), P4(0.35) e P2(0.35) onde f(x) = e x, e calcule o erro obtido pelas funções interpoladoras, de acordo com a estimativa de erro apresentada previamente. x 0 0.1 0.2 0.3 0.4 0.5 ex 1 1.11 1.22 1.35 1.49 1.65 Exercício: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p2(115) pela fórmula de Gregory-Newton: Interpolação de Gregory-Newton i 0 1 2 x 110 120 130 f(x) 2,041 2,079 2,114 18/02/2019 18 Exercício: Construir a tabela de diferenças ordinárias para a função f(x) tabelada abaixo e calcular p2(115) pela fórmula de Gregory-Newton: Resposta: Interpolação de Gregory-Newton i 0 1 2 x 110 120 130 f(x) 2,041 2,079 2,114 2 2 000015,000725,0425,1)( xxxp += 060,2)115(2 =p 06_Ajuste_de_Curvas_folheto.pdf 31/10/2018 1 Ajuste de Curvas Prof. Wellington Passos de Paula wpassos@gmail.com Programa 1. Introdução a) Caso Discreto 2. Método dos Quadrados Mínimos 3. Caso Não-Linear Ajuste de Curvas Introdução Introdução Anteriormente vimos uma forma de trabalhar com uma função definida por uma tabela: a interpolação polinomial Nem sempre a interpolação é aconselhável: Quando se quer aproximar um valor da função fora do intervalo de tabelamento (extrapolação). Quando os valores são medidas experimentais com erros. Nesse caso, a função deve passar pela barra de erros e não pelos pontos. Introdução Graficamente, a extrapolação e o ajuste por barras de erros são vistos abaixo: xexf )( x )(xf x xf Curva ajustada Curva extrapolada Barra de erros Introdução É necessário então ajustar essas funções tabeladas por uma função que seja uma “boa aproximação” e que nos permita “extrapolar” com certa margem de segurança Devemos então aproximar f(x) por outra função x), escolhida de uma família de funções ou por uma soma de funções em duas situações distintas: Caso discreto: quando a função f é dada por uma tabela de valores Caso contínuo: quando a função f é dada por sua forma analítica 31/10/2018 2 Introdução Caso discreto: Caso contínuo: Introdução – Caso Discreto Dados os pontos em um intervalo [a,b], devemos escolher funções , e constantes tais que a função se aproxime de Este modelo é dito linear pois os coeficientes a determinar aparecem linearmente Note que as funções podem ser funções não-lineares, como por exemplo )(,,......,)(,,)(, 2211 mm xfxxfxxfx )(,.......,)(,)( 21 xgxgxg n naaa ,.......,, 21 )(...)()()( 2211 xgxgxgx nnaaa )(xf naaa ,.......,, 21 )(,.......,)(,)( 21 xgxgxg n .......,1)(,)( 221 xxgexg x Introdução – Caso Discreto Problema 1 Como escolher as funções ? Podemos escolher as funções observando os pontos tabelados (diagrama de dispersão) ou a partir de conhecimentos teóricos do experimento Assim, dada uma tabela de pontos , devemos, primeiramente, colocar esses pontos em um gráfico cartesiano, o que vai nos permitir visualizar a curva que melhor se ajusta aos dados )(,...,)(,)( 21 xgxgxg n )(,,......,)(,,)(, 2211 mm xfxxfxxfx )(,.......,)(,)( 21 xgxgxg n Introdução – Caso Discreto Exemplo: Dada a tabela: Devemos construir o diagrama de dispersão x -1,0 -0,75 -0,6 -0,5 -0,3 0 0,2 0,4 0,5 0,7 1,0 f(x) 2,05 1,153 0,45 0,4 0,5 0 0,2 0,6 0,512 1,2 2,05 Introdução – Caso Discreto – Diagrama de Dispersão Introdução – Caso Discreto Escolhemos a partir da forma dos pontos no diagrama de dispersão Procuramos a função que se aproxime ao máximo de e que tenha a forma parábola passando pela origem Problema 2: Qual o valor de gera melhor ajuste da parábola? 2 1 )( xxg )(xf 2 11 )()( xxgx aa a 31/10/2018 3 Introdução – Caso Discreto Uma vez escolhidas as funções g1(x), g2(x), ..., gn(x), temos de estabelecer o conceito de proximidade entre as funções (x) e f(x) para obter as constantes a1, a2, a3, …, an Uma possibilidade é impor que o desvio entre f(x) e (x), ou seja, dk = (f(xk) - ϕ(xk)) seja mínimo para todos os pontos (k =1, 2, ...., m) Existem varias formas de impor que os desvios sejam mínimos. Estudaremos o Método dos Quadrados Mínimos Método dos Quadrados Mínimos Objetivo: encontrar os coeficientes aj tais que a função se aproxime ao máximo de f(x) Método dos Quadrados Mínimos - Consiste em escolher os aj’s de modo que a soma dos quadrados dos desvios seja mínima deve ser mínimo )()()()()( 2211 xgxxgxgx nnaaa m k kk m k k xxf 1 2 1 2 ))()((d Recapitulando, no caso discreto temos o tabelamento dos pontos como entrada do problema Dadas as funções , escolhidas de alguma forma, nosso objetivo então é encontrar os coeficientes tais que a função Se aproxime ao máximo de f(x) )(,...,,)(,,)(, 2211 mm xfxxfxxfx )(,.......,)(,)( 21 xgxgxg n )()()()( 2211 xgxgxgx nnaaa naaa ,.......,, 21 Método dos Quadrados Mínimos Seja o desvio em : Se a soma dos quadrados dos desvios é mínima, cada desvio será pequeno. Assim, aj’s devem ser tais que minimizem a função Se a aproximação (x) for perfeita somatório acima será nulo, que é o que acontece na interpolação )()( kkk xxf d m k kk m k k xxf 1 2 1 2 ))()((d )()( kkk xxf d m k kkn xxf 1 2 21 )]()([),,( aaa F kx Método dos Quadrados Mínimos Para obter um ponto mínimo devemos encontrar os números críticos, ou seja, aj’s tais que onde nj n j 2,1,0 ),,,( 21 aaa a F m k kkn xxf 1 2 21 )]()([),,,( aaa F Método dos Quadrados Mínimos m k knnkkk xgxgxgxf 1 2 2211 )]()()()([ aaa Calculando as derivadas, temos Igualando a zero, m k kjknnkkk j xgxgxgxgxf n 1 2211 ),,,( )]()][()()()([2 21 aaa a aaa F njxgxgxgxgxf m k kjknnkkk ,,2,1,0)]()][()()()([ 1 2211 aaa Método dos Quadrados Mínimos 31/10/2018 4 Ou seja, temos um sistema linear a resolver: 0)]()][()()()([ 0)]()][()()()([ 0)]()][()()()([ 1 2211 1 22211 1 12211 m k knknnkkk m k kknnkkk m k kknnkkk xgxgxgxgxf xgxgxgxgxf xgxgxgxgxf aaa aaa aaa Método dos Quadrados Mínimos Reescrevendo o sistema Sistema linear com n equações e com n incógnitas (a1, a2, a3, ..., an) m k knk m k nknkn m k knk m k kk m k nkkn m k kk m k kk m k nkkn m k kk xgxfxgxgxgxg xgxfxgxgxgxg xgxfxgxgxgxg 111 11 1 2 1 2 1 121 1 1 1 1 1 111 )()()]()([)]()([ )()()]()([)]()([ )()()]()([)]()([ aa aa aa Método dos Quadrados Mínimos O sistema linear pode ser reescrito na forma matricial Aa = b: onde A = (aij) tal que ou seja, A é uma matriz simétrica, e é tal que nnnnnn nn nn baaa baaa baaa aaa aaa aaa ... ... ... 2211 22222121 11212111 jiki m k kjkj m k kiij axgxgxgxga )()()()( 11 t n ]...,,,[ 21 aaaa t nbbbb ]...,,,[ 21 m k kiki xgxfb 1 )()( Método dos Quadrados Mínimos Método dos Quadrados Mínimos – Passo a Passo O funcionamento do Método dos Quadrados Mínimos pode ser dividido em 4 passos: Passo 1: Depois de escolhida a função ajuste (x) identificar nela as funções auxiliares g(x) tal que (x) seja do tipo: Passo 2: Montar o sistema de equações. O numero de equações do sistema é igual ao numero de funções auxiliares gi(x) ( igual ao numero de incógnitas ai ) )(...)()()( 2211 xgxgxgx nnaaa Método dos Quadrados Mínimos – Passo a Passo Passo 2: No caso da reta teremos um sistema com 2 equações: No caso de uma parábola teremos um sistema com 3 equações: exgxx 1)()( 121 aa xxg )(2 2 1 2 1 2221 1211 b b aa aa a a 2321)( xxx aaa 2 321 )()(,1)( xxgexxgxg 3 2 1 3 2 1 333231 232221 131211 b b b aaa aaa aaa a a a Método dos Quadrados Mínimos – Passo a Passo Passo 2: No caso de uma exponencial simples teremos um sistema com 1 equação: Passo 3: Calcular os coeficientes aij e bi do passo 2. Esses coeficientes são definidos pelos seguintes somatórios e após seu calculo obteremos números número de pontos experimentais xex 1)( a xexg )(1 1111 ba a jikj m k kiij axgxga )()( 1 m k kiki xgxfb 1 )()( 31/10/2018 5 Método dos Quadrados Mínimos – Passo a Passo Passo 4: Reescrever o sistema de equações do passo 2 (agora os aij e bi são números) e resolvê-lo, por exemplo, utilizando o método de eliminação de Gauss ou algum método iterativo (Gauss-Jacobi ou Gauss-Seidel). Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Nesse caso temos o que resulta em termos Para encontrarmos a1 e a2, resolveremos o sistema de 2 equações abaixo: . )(e1)( 21 xxgxg xxxf 21)()( aa . )(e1)( 21 xxgxg 2 1 2 1 2221 1211 b b aa aa a a Método dos Quadrados Mínimos - Retas Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Escrevendo o sistema em termos dos aij e bi, ficamos assim: . )(e1)( 21 xxgxg 4 1 22 4 1 221 4 1 12 4 1 12 4 1 211 4 1 11 )()()()()()( )()()()()()( k kk k kk k kk k kk k kk k kk xgxfxgxgxgxg xgxfxgxgxgxg aa aa Método dos Quadrados Mínimos - Retas Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Calculando os termos da primeira equação: = 1 = xk = 1 = 1 4 1 22222 1 4 1 11 41111)()()( k k k kk xgxgxg 913131211)()( 4 1 1 k kk xgxf 2218171512)()( 4 1 12 k kk xgxg . )(e1)( 21 xxgxg Método dos Quadrados Mínimos - Retas Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Calculando os termos da segunda equação: = 1 = xk = xk = xk = xk 5738372512)()( 4 1 2 k kk xgxf 1428752)()()( 2222 4 1 2 2 4 1 22 k k k kk xgxgxg 2281715121)()( 4 1 21 k kk xgxg . )(e1)( 21 xxgxg Método dos Quadrados Mínimos - Retas Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Logo, chegamos ao seguinte sistema de equações: Resolvendo o sistema acima pela Eliminação de Gauss chegamos à solução: . )(e1)( 21 xxgxg 57 9 14222 224 5714222 9224 2 1 21 21 a a aa aa T 14 5 7 2 a Método dos Quadrados Mínimos - Retas 31/10/2018 6 Exemplo: Encontre a reta de quadrados mínimos que melhor se ajusta aos pontos (2, 1), (5, 2), (7, 3), (8, 3). Calculemos para Solução: Substituindo os valores encontrados na equação original: . )(e1)( 21 xxgxg xxxf 21)()( aa xxxf 14 5 7 2 )()( Método dos Quadrados Mínimos - Retas Método dos Quadrados Mínimos - Retas Exemplo: Ajuste os dados abaixo pelo método dos quadrados mínimos utilizando uma parábola do tipo Solução: Para encontrarmos a1, a2 e a3 resolveremos o sistema de 3 equações a seguir: 2 321 2 321 )()(,1)()( xxgexxgxgxxx aaa x 1 2 3 4 5 6 7 8 f(x) 0,5 0,6 0,9 0,8 1,2 1,5 1,7 2,0 3 2 1 3 2 1 333231 232221 131211 b b b aaa aaa aaa a a a Método dos Quadrados Mínimos – Parábolas Exemplo: Ajuste os dados abaixo pelo método dos quadrados mínimos utilizando uma parábola do tipo Escrevendo o sistema em termos dos aij e bi, temos: 2 321 2 321 )()(,1)()( xxgexxgxgxxx aaa x 1 2 3 4 5 6 7 8 f(x) 0,5 0,6 0,9 0,8 1,2 1,5 1,7 2,0 8 1 33 8 1 332 8 1 231 8 1 13 8 1 23 8 1 322 8 1 221 8 1 12 8 1 13 8 1 312 8 1 211 8 1 11 )()()()()()()()( )()()()()()()()( )()()()()()()()( k kk k kk k kk k kk k kk k kk k kk k kk k kk k kk k kk k kk xgxfxgxgxgxgxgxg xgxfxgxgxgxgxgxg xgxfxgxgxgxgxgxg aaa aaa aaa Método dos Quadrados Mínimos – Parábolas 9,2 1 0,2 1 7,1 1 5,1 1 1,2 1 8,0 1 9,0 1 6,0 1 5,0)()( 8 1 1 k kk xgxf 361817 161514131211)()( 8 1 12 k kk xgxg 2041817 161514131211)()( 22 222222 8 1 13 k kk xgxg 8 1 222222222 1 8 1 11 811111111)()()( k k k kk xgxgxg Calculando os termos da primeira equação: = 1 = xk = 1 = xk 2 = 1 = 1 Método dos Quadrados Mínimos – Parábolas 368171 615141312111)()( 8 1 21 k kk xgxg 20487 654321)()()( 22 222222 8 1 2 2 8 1 22 k k k kk xgxgxg 50,5 8 0,2 7 7,1 6 5,1 5 1,2 4 8,0 3 9,0 2 6,0 1 5,0)()( 8 1 2 k kk xgxf 12968877 665544332211)()( 22 222222 8 1 23 k kk xgxg Calculando os termos da segunda equação: = 1 = xk = xk = xk = xk 2 = xk = xk Método dos Quadrados Mínimos – Parábolas 31/10/2018 7 2048171 615141312111)()( 22 222222 8 1 31 k kk xgxg 12968877 665544332211)()( 22 222222 8 1 32 k kk xgxg 319,1 8 0,2 7 7,1 6 5,1 5 1,2 4 8,0 3 9,0 2 6,0 1 5,0)()( 222 22222 8 1 3 k kk xgxf 877287 654321)()()( 44 444444 8 1 2 3 8 1 33 k k k kk xgxgxg Calculando os termos da terceira equação: = 1 = xk 2 = xk = xk 2 = xk 2 = xk 2 = xk 2 Método dos Quadrados Mínimos – Parábolas Logo, chegamos ao seguinte sistema de equações: Resolvendo o sistema acima pela Eliminação de Gauss chegamos à solução: 1,31987721296204 5,50129620436 2,9204368 321 321 321 aaa aaa aaa T0,03160,0871-0,7587a 1,319 5,50 2,9 87721296204 129620436 204368 3 2 1 a a a Método dos Quadrados Mínimos – Parábolas Exemplo: Ajuste os dados abaixo pelo método dos quadrados mínimos utilizando uma parábola do tipo Solução: Substituindo os valores encontrados na equação original: 2 321 2 321 )()(,1)()( xxgexxgxgxxx aaa x 1 2 3 4 5 6 7 8 f(x) 0,5 0,6 0,9 0,8 1,2 1,5 1,7 2,0 2 321)()( xxxxf aaa 20316.00871.07587,0)()( xxxxf Método dos Quadrados Mínimos – Parábolas Método dos Quadrados Mínimos – Parábolas Exemplo: Encontre a parábola através dos quadrados mínimos que melhor se ajusta aos pontos da tabela Solução: Nesse caso temos , ou seja, uma parábola. Todavia, pelo diagrama de dispersão (a seguir ) que uma parábola que passa pela origem seria uma boa escolha, logo o que resulta em termos Para encontrarmos a1, basta resolvermos a equação abaixo: x -1,0 -0,75 -0,6 -0,5 -0,3 0 0,2 0,4 0,5 0,7 1,0 f(x) 2,05 1,153 0,45 0,4 0,5 0 0,2 0,6 0,512 1,2 2,05 2 321)()( xxxxf aaa 2 1)()( xxxf a . )( 21 xxg 1111 ba a Método dos Quadrados Mínimos – Parábolas Introdução – Caso Discreto – Diagrama de Dispersão 31/10/2018 8 Montando os termos da equação: Calculando os termos: 11 1 1 11 1 111 )()()]()([ k kk k kk xgxfxgxg a 8464,2)1()7,0()5,0()4,0()2,0( )0()3,0()5,0()6,0(75,0)0,1()( 44444 444444 11 1 4 k kx Método dos Quadrados Mínimos – Parábolas 11 1 2 1 11 1 4 )()()( k kk k k xfxx a 21 11 1 11 11 1 2 1 )()()()( kk k kk k k xxgcomoxgxfxg a Calculando os termos: Substituindo na equação original, temos 8756,505,2)1(2,1)7,0( 512,0)5,0(6,0)4,0(2,0)2,0(0)0( 5,0)3,0(4,0)5,0(45,0)6,0( 153,175,005,2)0,1()()( 22 2222 222 22 11 1 2 k kk xfx 2 11 0642,2)(0642,28756,58464,2 xx aa Método dos Quadrados Mínimos – Parábolas Comentários: Note que a parábola pela origem que melhor ajusta os pontos fornecidos, através Método dos Quadrados Mínimos, é dada por Uma parábola da forma permite um melhor ajuste dos pontos, mas o sistema a ser resolvido é 3x3 com várias somas e produtos intermediários, o que aumenta o tempo de processamento 20642,2)( xx 2 321)( xxx aaa Método dos Quadrados Mínimos – Parábolas Exercício: Ajuste os dados abaixo pelo método dos quadrados mínimos utilizando uma parábola do tipo 2 321 2 321 )()(,1)()( xxgexxgxgxxx aaa x 0 1 2 3 4 f(x) 27 42 60 87 127 Método dos Quadrados Mínimos – Parábolas 0 0.5 1 1.5 2 2.5 3 3.5 4 20 40 60 80 100 120 140 Exercício: Ajuste os dados abaixo pelo método dos quadrados mínimos utilizando uma parábola do tipo Resposta: 2 321 2 321 )()(,1)()( xxgexxgxgxxx aaa x 0 1 2 3 4 f(x) 27 42 60 87 127 Método dos Quadrados Mínimos – Parábolas 221,464,702,28)( xxx 31/10/2018 9 Ajuste de Curvas Caso Não-Linear Prof. Wellington Passos de Paula wpassos@gmail.com Caso Não-Linear Em alguns casos, a família de funções escolhidas pode não ser linear nos parâmetros Ex: Função exponencial do tipo , sendo a1 e a2 positivos Para aplicar o Método dos Quadrados Mínimos torna-se necessário o uso de alguma transformação linear Ex: . Se e que é um problema linear nos parâmetros b1 e b2 xaexxf 21)()( a xyzey xa 211 )ln()ln( 2 aaa )ln( 11 ab )()ln( 2122 xxbbyb a Caso Não-Linear Aplicamos então o Método dos Quadrados Mínimos na resolução do problema linearizado. Utilizamos então os valores encontrados para calcular os parâmetros originais Observação: Os parâmetros a1 e a2 não serão ótimos dentro do critério dos quadrados mínimos, pois vamos aplicar o método ao problema linearizado e não ao problema original Caso Não-Linear Exemplo: Suponhamos que em um laboratório obtivemos experimentalmente os seguintes valores para f(x) sobre os pontos xi, i = 1, 2, ..., 8 Fazendo diagrama de dispersão dos dados x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 f(x) 36,547 17,264 8,155 3,852 1,820 0,860 0,406 0,246 Caso Não-Linear O gráfico de dispersão nos sugere um ajuste xexy 21)( aa Caso Não-Linear Como vimos, a lineaziração a ser feita é Logo, ajustaremos por quadrados mínimos, encontrando , onde . Assim, temos: )()ln()ln()ln( 211 2 xxeyz xa aaa )ln( yz xbbx 21)( eb )ln( 11 a 22 ab x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 31/10/2018 10 Solução: Montando o sistema linear. Utilizando (caso linear) e lembrando que b1 e b2 serão a solução desse sistema: 8 1 22 8 1 221 8 1 21 8 1 12 8 1 121 8 1 11 )()()()()()( )()()()()()( k kk k kk k kk k kk k kk k kk xgxfbxgxgbxgxg xgxfbxgxgbxgxg Caso Não-Linear x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 xxgxg )(e1)( 21 Solução: Calculando os termos da primeira equação: = 1 = xk = 1 = 1 8 1 2 1 8 1 11 8)()()( k k k kk xgxgxg 041,8)()( 8 1 1 k kk xgxz 3,0)()( 8 1 12 k kk xgxg Caso Não-Linear x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 Solução: Calculando os termos da segunda equação: = 1 = xk = xk = xk = xk 646,8)()( 8 1 2 k kk xgxz 59,3)()( 8 1 22 k kk xgxg 3,0)()( 8 1 21 k kk xgxg Caso Não-Linear x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 Solução: Logo, chegamos ao seguinte sistema de equações: Resolvendo o sistema acima pela Eliminação de Gauss chegamos à solução: 646.8 041,8 59,33,0 3,08 646,859,33,0 041,83,08 2 1 21 21 b b bb bb Tb 5,2099,1 Caso Não-Linear x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 Solução: Agora, Assim, a função Caso Não-Linear x -1,0 -0,7 -0,4 -0,1 0,2 0,5 0,8 1,0 z = ln(y) 3,599 2,849 2,099 1,349 0,599 -0,151 -0,901 -1,402 5,2 001,3)ln( 222 099,1 1111 1 aa aaa b eeb b xx eex 5,21 001,3)( 2 aa Existem outras situações em que será preciso fazer a linearização da curva ajustada aos pontos do diagrama de dispersão: Hipérbole: Curva Exponencial: Curva Trigonométrica: Caso Não-Linear )( 1 21 x x y aa x y z 21 1 aa )(21 xy x aa ))()ln()ln()ln(,0( 2121 xxbbxyzyse aa 1b 2b )()cos(21 xwxy aa ttwxt 21)()cos( aa 31/10/2018 11 Exercício: O número de bactérias, por unidade de volume, existente em uma cultura após x horas é dado na tabela abaixo: a) Ajuste os dados acima à curva pelo método dos quadrados mínimos. Resposta: b) Quantas horas seriam necessárias para que o número de bactérias por unidade de volume ultrapasse 2000? Resposta: 11,64 hrs Caso Não-Linear x 0 1 2 3 4 5 6 y 32 47 65 92 132 190 275 xey 21 aa xeyxy 355,0104,32355,0469,3ln 07_Integracao_folheto.pdf 16/05/2019 1 Integração Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Programa 1. Introdução 2. Fórmulas de Newton – Cotes a) Regra dos Trapézios b) Regra 1/3 de Simpson c) Regra 3/8 de Simpson 3. Quadratura Gaussiana Integração Introdução Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Introdução O Cálculo Diferencial e Integral ensina que se y = f(x) é uma função contínua em [a, b] , então para se obter basta determinar uma primitiva, isto é, uma função F(x), tal que F’(x) = f(x), de forma que b a dxf(x) I F(a)F(b)dxf(x) I b a Introdução Problemas Pode não ser fácil, ou impossível, expressar F(x) por meio de uma combinação finita de funções elementares Há situações nas quais y = f(x) é conhecida apenas em um conjunto discreto de pontos Solução Nessas situações, avalia-se numericamente! b a dxf(x)I Introdução Estratégia Para calcular numericamente vamos expressar f(x) como um polinômio no intervalo [a,b]. Deduziremos expressões que têm a forma onde Quando escrevemos uma integral na forma acima , estamos implementando o formalismo de Newton - Cotes dxxf b a )( )(....)()()( 1100 nn b a xfAxfAxfAdxxf ....,,2,1,0 com , nibaxi 16/05/2019 2 Integração Fórmulas de Newton - Cotes Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Fórmulas de Newton - Cotes No procedimento de Newton - Cotes o polinômio aproxima f(x) em pontos de [a,b], igualmente espaçados. Se os subintervalos têm comprimento h, então as fórmulas fechadas de Newton - Cotes para integração têm a forma n ab hxx xfAxfAxfAxfA dxxfdxxf ii ii n i nn x x b a n )( onde )(....)()( )()( 1 0 1100 0 Fórmulas de Newton - Cotes Os coeficientes Ai das formulas fechadas de Newton - Cotes são determinados de acordo com o grau do polinômio aproximador de f(x) As fórmulas abertas de Newton - Cotes, construídas de forma análoga às fechadas, diferem pelo fato que Vamos estudar as seguintes fórmulas fechadas de Newton – Cotes: Regra dos Trapézios Regra 1/3 de Simpson Regra 3/8 de Simpson baxx n , e 0 Integração Regra dos Trapézios Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Regra dos Trapézios A estratégia da Regra dos Trapézios é aproximar a função f(x) por um polinômio de ordem 1 (reta) Veremos que, nessa aproximação a integral da função f(x) pode ser aproximada pela área de um trapézio ,, e 10 baxx Regra dos Trapézios Utilizando a Fórmula de Lagrange para expressar p1(x), que interpola f(x) em teremos o seguinte: sendo Fazendo h = (x1 – x0)/n, onde nesse caso n = 1 e substituindo os fatores de Lagrange no polinômio podemos reescrevê-lo assim: ,, e 10 baxx )()()()()( 11001 xLxfxLxfxp 01 0 1 10 1 0 )()( xx xx xLe xx xx xL h xx xf h xx xfxp 01 1 01 )()()( 16/05/2019 3 T x x Idxxf h xx xf h xx )()( 1 0 0 11 0 )()( 2 10 xfxf h IT Regra dos Trapézios Pela nossa aproximação, a integral da função f(x) será escrita por: dxxpdxxf xb xa b a )()( 1 1 0 Note que IT é a área do trapézio de altura h = x1 - x0 e de base f(x0) e f(x1) Regra dos Trapézios )x(fy 0y 1xb 0xa )x(py 1y Regra dos Trapézios Ao substituir a área sob a curva f(x) pela área do trapézio estamos realizando uma aproximação e cometendo um erro. Verifica-se que este erro é dado por: ],[ onde )(max 12 3 baccf ab ET Regra dos Trapézios Exemplo: Considere a integral a) Calcule a integral utilizando a Regra dos Trapézios. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Solução: a) Utilizando a equação de IT, , e, sendo a = 1, b = 7, h = 6, temos: dx x2 7 1 1 )()( 2 10 xfxf h IT 22 7 1 1 1 2 6 TI 0612245,3 Regra dos Trapézios Exemplo: Considere a integral a) Calcule a integral utilizando a Regra dos Trapézios. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Solução: b) Para o erro, temos que: Como a derivada segunda de f(x) é f’’(x) = 6x-4, logo: dx x2 7 1 1 ],[ )(max 12 6 3 baccfET 6 12 63 108 TE )1(''12 6 3 fET Regra dos Trapézios Exercício: Considere a integral a) Calcule a integral utilizando a Regra dos Trapézios. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 1,8591 b) dxex 1 0 0,2265TE 16/05/2019 4 Regra dos Trapézios Exercício: Considere a integral a) Calcule a integral utilizando a Regra dos Trapézios. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 32 b) dxx 56 9 1 384TE Regra dos Trapézios Repetida A Regra dos Trapézios é uma aproximação um pouco grosseira para o valor da integral (erros muito grandes) Afim de melhorar essa aproximação podemos fazer várias subdivisões no intervalo [a, b] e aplicar a Regra dos Trapézios repetidas vezes, conforme a figura abaixo: Regra dos Trapézios Repetida Dividindo o intervalo [a, b] em n subdivisões iguais de largura h, sendo h calculado da seguinte forma: Logo, os valores de cada ponto xi das subdivisões pode se obtidos através da expressão: Podemos então escrever a integral de f(x) como sendo a soma das áreas dos n trapézios pequenos contidos dentro do intervalo [a, b]: tal que Ai = área do trapézio i, com i = 1, 2, …, n. n ab h )( hixxi 0 n b a AAAAdxxf ...)( 321 Regra dos Trapézios Repetida Sendo Temos então que o valor numérico da integral calculada segundo a Regra dos Trapézios Repetida será ou, na forma resumida: )()( 2 1 iii xfxf h A )()( 2 )()( 2 ...)()( 2 )()( 2 )()( 2 )( 11232 2110 nnnn b a xfxf h xfxf h xfxf h xfxf h xfxf h dxxf TR n i in b a Ixfxfxf h dxxf 1 1 0 )(2)()( 2 )( Regra dos Trapézios Repetida A estimativa para o erro gerado pela aplicação da Regra dos Trapézios Repetida é dada por: sendo n o número de subdivisões no intervalo [a, b]. Comparando com o erro na Regra dos Trapézios Simples, temos: ou seja: ],[ onde )(max 12 )( 2 3 baccf n ab ETR )(max 12 )( )(max 12 2 33 cf n ab Ecf ab E TRT 2n E E TTR Regra dos Trapézios Repetida Se quisermos fazer o caminho inverso, ou seja, calcular quantas subdivisões são necessárias para atingir um certo valor de erro, devemos fazer o seguinte cálculo: )(max 12 )( n 3 cf E ab TR 16/05/2019 5 Regra dos Trapézios Repetida Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando 10 subintervalos e a Regra dos Trapézios Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-3? dx x2 7 1 1 Regra dos Trapézios Repetida Solução: a) Fazendo 10 subintervalos no intervalo [1, 7] temos h = 0,6 Como , logo: x0 = 1, x1 = 1,6 , x2 = 2,2 , x3 = 2,8 , x4 = 3,4 , x5 = 4,0 , x6 = 4,6 , x7 = 5,2 , x8 = 5,8 , x9 = 6,4, x10 = 7. Aplicando então a Regra dos Trapézios Repetida: hixxi 0 1 1 0 )(2)()( 2 n i inTR xfxfxf h I 2 9 2 8 2 7 2 6 2 5 2 4 2 3 2 2 2 1 2 10 2 0 111111111 2 11 2 xxxxxxxxx xx h ITR Regra dos Trapézios Repetida Solução: a) Fazendo 10 subintervalos no intervalo [1, 7] temos h = 0,6 Como , logo: x0 = 1, x1 = 1,6 , x2 = 2,2 , x3 = 2,8 , x4 = 3,4 , x5 = 4,0 , x6 = 4,6 , x7 = 5,2 , x8 = 5,8 , x9 = 6,4, x10 = 7. Aplicando então a Regra dos Trapézios Repetida: hixxi 0 222222222 22 4,6 1 8,5 1 2,5 1 6,4 1 0,4 1 4,3 1 8,2 1 2,2 1 6,1 1 2 7 1 1 1 2 6,0 TRI 0,9134TRI Regra dos Trapézios Repetida Solução: b) Para o erro, temos: Como a derivada segunda de f(x) é f’’(x) = 6x-4, logo: )(max 12 )( 2 3 cf n ab ETR 6 1012 1) -(7 2 3 08,1 TRE )(max 12 )( 2 3 cf n ab ETR Regra dos Trapézios Repetida Solução: c) O número de subdivisões para que o erro fosse menor do que 10-3 pode ser obtido por: 6 1012 )17( 3 3 329n )(max 12 )( n 3 cf E ab TR 63,328n Regra dos Trapézios Repetida Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando 10 subintervalos e a Regra dos Trapézios Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-3? Respostas: a) 1,7197 b) c) 15,051 16 subdivisões dxex 1 0 0,00227TRE 16/05/2019 6 Regra dos Trapézios Repetida Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando 8 subintervalos e a Regra dos Trapézios Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-4? Respostas: a) 37,818 b) c) 3,84 x 106 subdivisões dxx 56 9 1 6TRE Integração Regra 1/3 de Simpson Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Regra 1/3 de Simpson Consideremos agora que se queira aproximar f(x) por um polinômio interpolador de ordem 2 (parábola) Regra 1/3 de Simpson Utilizando a Fórmula de Lagrange para expressar p2(x): com , e temos ainda que: , assim: x0 – x1 = -h x0 – x2 = -2h x1 – x0 = h x1 – x2 = -h x2 – x0 = 2h x2 – x1 = h )()()()()()()( 2211002 xfxLxfxLxfxLxp 2010 21 0 )( xxxx xxxx xL 2101 20 1 )( xxxx xxxx xL 1202 10 2 )( xxxx xxxx xL bhxx 202e e a 010 hxxx Regra 1/3 de Simpson Logo, S bx ax b a Idxxpdxxf )()( 2 2 0 )( 2 )()( 2 )( 2 10 1 20 0 21 2 xf hh xxxx xf hh xxxx xf hh xxxx xp dxxxxx h xf dxxxxx h xf dxxxxx h xf I x x x x x x S 102 2 202 1 212 0 2 0 2 0 2 0 2 )( )( 2 )( )()(4)( 3 210 xfxfxf h IS Regra 1/3 de Simpson De modo análogo à Regra dos Trapézios, na Regra 1/3 de Simpson estamos realizando uma aproximação e cometendo um erro. Verifica-se que este erro é dado por como h = (b – a) / 2 h5 = (b – a)5 / 32, tem-se: ],[ onde )(max 90 4 5 baccf h ES ],[ onde )(max 2880 )( 4 5 baccf ab ES 16/05/2019 7 Regra 1/3 de Simpson Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 1/3 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. dx x2 7 1 1 Regra 1/3 de Simpson Solução: a) Como agora temos n = 2 subdivisões dentro do intervalo [a, b], teremos h = (b – a)/2 = (7 – 1)/2 = 3. Logo: x0 = 1, x1 = x0 + h = 4, x2 = 7. Assim o valor da integral será: 1,27 )()(4)( 3 210 xfxfxf h IS 222 7 1 4 1 4 1 1 3 3 SI Regra 1/3 de Simpson Solução: b) Para o erro, temos: Como a derivada quarta de f(x) é f4(x) = 120x-6 logo )(max 2880 )( 4 5 cf ab ES 324 TRE )(max 2880 )( 4 5 cf ab ES 120 2880 6 5 Regra 1/3 de Simpson Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 1/3 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 1,719 b) dxex 1 0 410439,9 SE Regra 1/3 de Simpson Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 1/3 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 37,3333 b) dxx 56 9 1 13824SE Regra 1/3 de Simpson Repetida Vamos agora repetir o procedimento anterior para n pares de subintervalos. Definimos o número de subintervalos pela letra m = 2n Por que precisamos definir a variável n? A cada par de subintervalos temos 3 pontos para ajustar uma parábola ( p2(x) ) 16/05/2019 8 Regra 1/3 de Simpson Repetida Tomemos para m = 2n m é par Aplica-se a Regra 1/3 de Simpson repetidas vezes no intervalo [a, b] = [x0, xm]: tal que Ai = área do subintervalo i, com i = 1, 2, …, n. mihxx m ab h ii ...,,2,1,0 com 1 n b a AAAAdxxf ...)( 321 Regra 1/3 de Simpson Repetida Sendo Temos então que o valor numérico da integral calculada segundo a Regra 1/3 de Simpson Repetida será )()(4)( 3 2101 xfxfxf h A )()(4)( 3 ...)()(4)( 3 )()(4)( 3 )()( 12432 210 0 mmm x x b a xfxfxf h xfxfxf h xfxfxf h dxxfdxxf m ))](...)()((4 ))(...)()((2)()([ 3 )( 131 2420 m mm b a xfxfxf xfxfxfxfxf h dxxf Regra 1/3 de Simpson Repetida Logo: onde: Valor da função nos pontos de índices pares dentro do intervalor [a, b], excluindo as extremidades Valor da função nos pontos de índices ímpares dentro do intervalor [a, b], excluindo as extremidades SR m i i m i im b a Ixfxfxfxf h dxxf 2 1 12 1 2 1 20 )(4)(2)()( 3 )( 1 2 1 2 )(2 m i ixf 2 1 12 )(4 m i ixf Regra 1/3 de Simpson Repetida A estimativa para o erro gerado pela aplicação da Regra 1/3 de Simpson Repetida é dada por: como h = (b – a) / m h5 = (b – a)5 / 32n5, tem-se: ],[ onde )(max 90 4 5 baccf h nESR ],[ onde )(max 2880 )( 4 4 5 baccf n ab ESR Regra 1/3 de Simpson Repetida Comparando com o erro na Regra 1/3 de Simpson Simples, temos: ou seja: Se quisermos fazer o caminho inverso, ou seja, calcular quantas subdivisões são necessárias para atingir um certo valor de erro, devemos fazer o seguinte cálculo: 4 4 5 )(max 2880 )( n cf E ab SR )(max 2880 )( )(max 2880 )( 4 4 5 4 5 cf n ab Ecf ab E SRS 4n E E SSR Regra 1/3 de Simpson Repetida Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando 10 subintervalos e a Regra 1/3 de Simpson Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-3? dx x2 7 1 1 16/05/2019 9 Regra 1/3 de Simpson Repetida Solução: a) Fazendo 10 subintervalos no intervalo [1, 7] temos h = 0,6. Como , logo: x0 = 1 , x1 = 1,6 , x2 = 2,2 , x3 = 2,8 , x4 = 3,4 , x5 = 4 , x6 = 4,6 , x7 = 5,2 , x8 = 5,8 , x9 = 6,4 , x10 = 7. Aplicando então a Regra 1/3 de Simpson Repetida: hixxi 0 2 1 12 1 2 1 20 )(4)(2)()( 3 m i i m i imSR xfxfxfxf h I Regra 1/3 de Simpson Repetida Solução: a) Calculando os somatórios: )()()()()()( 8642 41 2 10 1 2 1 2 1 2 xfxfxfxfxfxf i i m i i )()()()()()()( 97531 5 2 10 1 12 2 1 12 xfxfxfxfxfxfxf i i m i i 3701,0 8,5 1 6,4 1 4,3 1 2,2 1 2222 642,0 4,6 1 2,5 1 4 1 8,2 1 6,1 1 22222 Regra 1/3 de Simpson Repetida Solução: a) Logo: 2 1 12 1 2 1 2100 )(4)(2)()( 3 m i i m i iSR xfxfxfxf h I 8657,06427,04701,02 7 1 1 1 3 6,0 22 SRI Regra 1/3 de Simpson Repetida Solução: b) Para o erro, temos: Como a derivada quarta de f(x) é f4(x) = 120x-6 logo )(max 2880 )( 4 4 5 cf n ab ESR 5184,0 SRE )(max 2880 )( 4 4 5 cf n ab ESR 120 52880 6 4 5 Regra 1/3 de Simpson Repetida Solução: c) O número de subdivisões para que o erro fosse menor do que 10-3 pode ser obtido por: Logo, precisamos de 52 subdivisões! 4 3 5 120 102880 )17( 4 4 5 )(max 2880 )( n cf E ab SR 85,25n 26n Regra 1/3 de Simpson Repetida Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando 8 subintervalos e a Regra 1/3 de Simpson Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-10? Respostas: a) 1,7183 b) c) n = 56 m = 112 subdivisões dxex 1 0 -6103,6869SRE 16/05/2019 10 Integração Regra 3/8 de Simpson Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Regra 3/8 de Simpson Consideremos agora que se queira aproximar f(x) por um polinômio interpolador de ordem 3 Utilizando a Fórmula de Lagrange para expressar p3(x) que interpola f(x) nos pontos: bhxxhxxhxxx 3,2,,a 0302010 Regra 3/8 de Simpson Utilizando a Fórmula de Lagrange para expressar p3(x): com e )()()()()()()()()( 332211003 xfxLxfxLxfxLxfxLxp ,)( 302010 321 0 xxxxxx xxxxxx xL ,)( 312101 320 1 xxxxxx xxxxxx xL 321202 310 2 )( xxxxxx xxxxxx xL 231303 210 3 )( xxxxxx xxxxxx xL Regra 3/8 de Simpson Logo, 8/33 )()( 3 0 S bx ax b a Idxxpdxxf )( 23 )( 2 )( 2 )( 32 )( 3 210 2 310 1 320 0 321 3 xf hhh xxxxxx xf hhh xxxxxx xf hhh xxxxxx xf hhh xxxxxx xp Regra 3/8 de Simpson Logo, dxxxxxxx h xf dxxxxxxx h xf dxxxxxxx h xf dxxxxxxx h xf I x x x x x x x xS 2103 3 3103 2 3203 1 3213 0 8/3 3 0 3 0 3 0 3 0 6 )( 2 )( 2 )( 6 )( )()(3)(3)( 8 3 32108/3 xfxfxfxfhIS Regra 3/8 de Simpson O erro cometido na aproximação pela Regra 3/8 de Simpson é dado por: como h = (b – a) / 3 h5 = (b – a)5 / 243, tem-se: ],[ onde )(max 80 3 45 8/3 baccfhES ],[ onde )(max 19440 )( 4 5 8/3 baccf ab ES 16/05/2019 11 Regra 3/8 de Simpson Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 3/8 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. dx x2 7 1 1 22228/3 7 1 5 1 3 3 1 3 1 1 2 8 3 SI Regra 3/8 de Simpson Solução: a) Como agora temos n = 3 subdivisões dentro do intervalo [a, b], teremos h = (b – a)/3 = (7 – 1)/3 = 2. Logo: x0 = 1, x1 = x0 + h = 3, x2 = x1 + h = 5, x3 = 7. Assim o valor da integral será: 1,105 )()(3)(3)( 8 3 32108/3 xfxfxfxfhIS Regra 3/8 de Simpson Solução: b) Para o erro, temos: Como a derivada quarta de f(x) é f4(x) = 120x-6 logo )(max 19440 )( 4 5 8/3 cf ab ES 48 TRE )(max 19440 )( 4 5 8/3 cf ab ES 120 19440 6 5 Regra 3/8 de Simpson Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 3/8 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 1,716 b) dxex 1 0 410398,1 SE Regra 3/8 de Simpson Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Regra 3/8 de Simpson. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. Respostas: a) 37,610 b) dxx 56 9 1 0482SE Regra 3/8 de Simpson Repetida Vamos agora repetir o procedimento anterior para n trios de subintervalos. Definimos o número de subintervalos pela letra m = 3n Por que precisamos definir a variável n? A cada trio de subintervalos temos 4 pontos para ajustar uma função de grau 3 ( p3(x) ) 16/05/2019 12 Regra 3/8 de Simpson Repetida Tomemos para m = 3n m é múltiplo de 3 Aplica-se a Regra 3/8 de Simpson repetidas vezes no intervalo [a, b] = [x0, xm]: tal que Ai = área do subintervalo i, com i = 1, 2, …, n. mihxx m ab h ii ...,,2,1 com 1 n b a AAAAdxxf ...)( 321 Regra 3/8 de Simpson Repetida Sendo Temos então que o valor numérico da integral calculada segundo a Regra 3/8 de Simpson Repetida será )()(3)(3)( 8 3 32101 xfxfxfxfhA )()(3)(3)( 8 3 )()(3)(3)( 8 3 ...)()(3)(3)( 8 3 )()(3)(3)( 8 3 )()( 123 3456 6543 3210 0 mmmm mmmm x x b a xfxfxfxfh xfxfxfxfh xfxfxfxfh xfxfxfxfhdxxfdxxf m Regra 3/8 de Simpson Repetida Temos então que o valor numérico da integral calculada segundo a Regra 3/8 de Simpson Repetida será Logo: ))]()(...)()()()((3 ))(...)()((2)()([ 8 3 )( 125421 3630 mm mm b a xfxfxfxfxfxf xfxfxfxfxfhdxxf 3 1 13 3 1 23 1 3 1 308/3 )(3)(3)(2)()( 8 3 m i i m i i m i imSR xfxfxfxfxfhI Regra 3/8 de Simpson Repetida A estimativa para o erro gerado pela aplicação da Regra 3/8 de Simpson Repetida é dada por: como h = (b – a) / m h5 = (b – a)5 / 243n5, tem-se: ],[ onde )(max 80 3 45 8/3 baccfhnESR ],[ onde )(max 19440 )( 4 4 5 8/3 baccf n ab ESR Regra 3/8 de Simpson Repetida Comparando com o erro na Regra 3/8 de Simpson Simples, temos: ou seja: Se quisermos fazer o caminho inverso, ou seja, calcular quantas subdivisões são necessárias para atingir um certo valor de erro, devemos fazer o seguinte cálculo: 4 4 5 )(max 19440 )( n cf E ab SR )(max 19440 )( )(max 19440 )( 4 4 5 8/3 4 5 8/3 cf n ab Ecf ab E SRS 4 8/3 8/3 n E E SSR Regra 3/8 de Simpson Repetida Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando 12 subintervalos e a Regra 3/8 de Simpson Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-3? dx x2 7 1 1 16/05/2019 13 Regra 3/8 de Simpson Repetida Solução: a) Fazendo 12 subintervalos no intervalo [1, 7] temos h = 0,5 Como , logo: x0 = 1 , x1 = 1,5 , x2 = 2,0 , x3 = 2,5 , x4 = 3,0 , x5 = 3,5 , x6 = 4,0 , x7 = 4,5 , x8 = 5,0 , x9 = 5,5 , x10 = 6 , x11 = 6,5 , x12 = 7. Aplicando então a Regra 3/8 de Simpson Repetida: hixxi 0 3 1 13 3 1 23 1 3 1 308/3 )(3)(3)(2)()( 8 3 m i i m i i m i imSR xfxfxfxfxfhI Regra 3/8 de Simpson Repetida Solução: a) Calculando os somatórios: )()()()()(2 963 31 3 12 1 3 1 3 1 3 xfxfxfxfxf i i m i i )()()()()()(3 10741 4 3 12 1 23 3 1 23 xfxfxfxfxfxf i i m i i 256,0 5,5 1 4 1 5,2 1 222 633,0 6 1 5,4 1 3 1 5,1 1 2222 )()()()()()(3 11852 4 3 12 1 13 3 1 13 xfxfxfxfxfxf i i m i i 395,0 5,6 1 5 1 5,3 1 2 1 2222 Regra 3/8 de Simpson Repetida Solução: a) Logo: 395,03633,03256,02 7 1 1 1 5,0 8 3 228/3SR I 3 1 13 3 1 23 1 3 1 31208/3 )(3)(3)(2)()( 8 3 m i i m i i m i iSR xfxfxfxfxfhI 889,0 Regra 3/8 de Simpson Repetida Solução: b) Para o erro, temos: Como a derivada quarta de f(x) é f4(x) = 120x-6 logo )(max 19440 )( 4 4 5 8/3 cf n ab ESR 0,1875 8/3 SRE )(max 19440 )( 4 4 5 8/3 cf n ab ESR 120 419440 6 4 5 Regra 3/8 de Simpson Repetida Solução: c) O número de subdivisões para que o erro fosse menor do que 10-3 pode ser obtido por: Logo, precisamos de 45 subdivisões! 4 3 5 120 1019440 )17( 4 4 5 )(max 19440 )( n cf E ab SR 80,14n 51n Regra 3/8 de Simpson Repetida Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando 12 subintervalos e a Regra 3/8 de Simpson Repetida. b) Calcule a estimativa para o erro gerado na utilização dessa técnica numérica. c) Qual é o número mínimo de subdivisões, de modo que o erro seja inferior a 10-10? Respostas: a) 1,715 b) c) n = 35 m = 105 subdivisões! dxex 1 0 -7105,462SRE 16/05/2019 14 Comentários Percebemos que à medida que aumentamos o grau do polinômio, a convergência do método é mais fácil As demais fórmulas fechadas de integração de Newton - Cotes trabalham com polinômios de graus n=4, n=5, ... Para um n qualquer, a fórmula de Newton – Cotes é apresentada no próximo slide Comparativo Fórmula de Newton - Cotes para n qualquer )( 0 nx x ii dxxLA dxxLxfxLxfxLxf nn x x n )()(...)()()()( 1100 0 )()(...)()()()( 000 1100 nnn x x nn x x x x dxxLxfdxxLxfdxxLxf Lagrange de forma utilizando )()( 00 dxxpdxxf n x x xb xa nn )(....)()( 1100 nn xfAxfAxfA Integração Quadratura Gaussiana Prof. Wellington Passos de Paula wpassos@ufsj.edu.br Quadratura Gaussiana Os métodos de integração numérica apresentados até agora (a saber, Regra dos Trapézios e de Simpson) tomam como base uma regra simples para escolha dos pontos de avaliação da função f(x), onde Esses métodos são particularmente adequados para dados tabulados de forma regular, tais como algumas medidas de laboratório e valores obtidos de programas de computador que produzem tabelas hxx ii 1 Quadratura Gaussiana Se, por outro lado, tivermos a liberdade de escolher os pontos nos quais a função f(x) é avaliada, então uma escolha cuidadosa pode levar a uma maior precisão da avaliação da função Essa é a proposta do método conhecido por Quadratura Gaussiana, Integração Gaussiana ou Integração de Gauss-Legendre Quadratura Gaussiana Outra vantagem da Quadratura Gaussiana no erro relativo à integração: Os erros dos poliômios interpoladores gerados pelas fórmulas de Newton-Cotes envolvem a (n+1)-ésima ou (n+2)-ésima derivadas. Logo, essas integrações são exatas para polinômios de grau < n+1 ou < n+2, respectivamente. A Quadratura de Gauss integra exatamente polinômios de grau < 2n+2 16/05/2019 15 Quadratura Gaussiana Partimos então da fórmula para a integral utilizada nos métodos de Newton - Cotes onde os coeficientes Ai e os pontos xi para i = 0, 1, 2, .., n devem ser determinados de modo a obter a melhor precisão possível Diferença Importante: Uso de Partições Não-Regulares nn b a xfAxfAxfAdxxfI ....)( 1100 Quadratura Gaussiana Comecemos então o desenvolvimento para dois pontos: Por simplicidade tomemos o intervalo [-1, 1]. Note que sempre é possível passar do intervalo [a, b] --> [-1, 1] através da transformação: 1100)( xfAxfAdxxfI b a dtabdttxdx )( 2 1 )( 1,1 para )( 2 1 )( 2 1 )( tabtabtx Segue devemos então encontrar valores de A0, A1, t0, t1 que tornem a exatidão da integral a maior possível. Para isto o método é construído para ser exato para polinômios de graus até 3, pois deste modo teremos quatro incógnitas (A0, A1, t0, t1) e quatro equações )( 2 1 )( 2 1 )( 2 1 )( onde )( 1 1 ababtabftF dttF Quadratura Gaussiana dxxfI b a )( dttxtxf )())(( 1 1 1100 1 1 )( tFAtFAdttFI Quadratura Gaussiana Sejam Note que qualquer polinômio de grau 3 é combinação das funções acima. Assim, impomos que a fórmula da Quadratura Gaussiana seja exata para estes polinômios, segue: 3 3 2 210 )(,)(,)(,1)( ttFttFttFtF )()()()()( 332211003 tFatFatFatFatP )()()()( )()()()( 13103031210202 11101011010000 tFAtFAatFAtFAa tFAtFAatFAtFAa dttFadttFa dttFadttFa 1 1 33 1 1 22 1 1 11 1 1 00 )()( )()( dttP )(3 1 1 Quadratura Gaussiana )()( 131030 tPAtPA )()()()( )()()()( 1313121211111010 0303020201010000 tFAatFAatFAatFAa tFAatFAatFAatFAa ))()()()(( ))()()()(( 1331221111001 0330220110000 tFatFatFatFaA tFatFatFatFaA Exatamente a fórmula da Quadratura Gaussiana quando temos somente 2 pontos Quadratura Gaussiana Considerando podemos determinar as incógnitas A0, A1, t0, t1 através de ou seja, igualamos a equação à integral analítica de , gerando então um sistema não- linear 4 x 4. Vejamos: 3 3 2 210 )(,)(,)(,1)( ttFttFttFtF 3,2,1,0 para )( 1100 1 1 1 1 ktAtAdttdttFI kkk k kk tAtA 1100 )(tFk 16/05/2019 16 Quadratura Gaussiana Obtemos o sistema 20 10 0 11 0 00 1 1 0 AAtAtAdttk 01 1100 1 11 1 00 1 1 1 tAtAtAtAdttk 3/22 2 11 2 00 2 01 2 00 1 1 2 tAtAtAtAdttk 03 3 11 3 00 3 11 3 00 1 1 3 tAtAtAtAdttk Quadratura Gaussiana Resolvendo o sistema, obtemos de modo que podemos escrever a Fórmula de Quadratura Gaussiana, que é exata para polinômios de graus até 3, como 3/3e1 1010 ttAA 3 3 3 3 )( 1 1 FFdttFIGauss Quadratura Gaussiana Para 3 pontos, a fórmula da quadratura gaussiana é exata para polinômios de graus até 5. Então, Analogamente, qualquer polinômio de grau 5 pode ser escrito em termos de 221100 1 1 )( tFAtFAtFAdttFI 5 5 4 4 3 3 2 210 )(,)(,)(,)(,)(,1)( ttFttFttFttFttFtF Quadratura Gaussiana Agora podemos determinar as incógnitas A0, A1, A2, t0, t1, t2 através do sistema linear 6 x 6 abaixo: Escrevendo explicitamente o sistema: 5,4,3,2,1,0 para )( 221100 1 1 1 1 k tAtAtAdttdttFI kkkk k Quadratura Gaussiana Obtemos: 20 210 0 22 0 11 0 00 1 1 0 AAAtAtAtAdttk 01 221100 1 22 1 11 1 00 1 1 1 tAtAtAtAtAtAdttk 3/22 2 21 2 11 2 00 2 22 2 11 2 00 1 1 2 tAtAtAtAtAtAdttk 05 5 22 5 11 5 00 5 22 5 11 5 00 1 1 5 tAtAtAtAtAtAdttk 03 3 22 3 11 3 00 3 22 3 11 3 00 1 1 3 tAtAtAtAtAtAdttk 5/24 4 21 4 11 4 00 4 22 4 11 4 00 1 1 4 tAtAtAtAtAtAdttk Quadratura Gaussiana Resolvendo o sistema, obtemos de modo que podemos escrever a Fórmula de Quadratura Gaussiana, que é exata para polinômios de graus inferiores a 5, como 0e 5 3 , 9 8 e 9 5 120120 tttAAA 5 3 9 5 0 9 8 5 3 9 5 )( 1 1 FFFdttFIGauss 16/05/2019 17 Quadratura Gaussiana A estimativa para o erro gerado pela aplicação da Quadratura Gaussiana é dada por: onde n é o grau do polinômio interpolado ],[ onde )(max ])!22)[(32( ])!1[(2 )22( 3 432 baccf nn n E n n QG Quadratura Gaussiana Exemplo: Calcule utilizando o Método da Quadratura Gaussiana para 2 e 3 pontos. Solução: Temos no intervalo [1, 3]. Fazendo a mudança de variáveis dxeI x3 3 1 xexf 3)( 1,1 temos3,1 para 2)( 2 1 )( 2 1 )( tx tabtabtx 2t3e dtdtabdttxdx 1)( 2 1 )( )( 2 1 )( 2 1 )( 2 1 )( ababtabftF Quadratura Gaussiana Exemplo: Calcule utilizando o Método da Quadratura Gaussiana para 2 e 3 pontos. Solução: Para 2 pontos: dxeI x3 3 1 1018,523 :Exato 3 1 dxeI x 3 3 3 3 )(3 1 1 3 1 FFdttFIdxe Gauss x 2 3 3 2 3 3 33 ee 9309,51 Quadratura Gaussiana Exemplo: Calcule utilizando o Método da Quadratura Gaussiana para 2 e 3 pontos. Solução: Para 3 pontos: dxeI x3 3 1 1018,523 :Exato 3 1 dxeI x dttFIdxe Gauss x )(3 1 1 3 1 25 3 2 25 3 9 5 3 9 8 3 9 5 3 eee 5 3 9 5 0 9 8 5 3 9 5 FFF 1004,52 Quadratura Gaussiana A tabela abaixo compara o Método da Quadratura Gaussiana com a Regra 1/3 de Simpson para dxeI x3 3 1 Exato Gauss 2 pontos Gauss 3 pontos Simpson 3 pontos Simpson 5 pontos Simpson 7 pontos Valor 52,1018 51,9309 52,1004 52,3601 52,1194 52,1053 Erro 0,1709 0,0014 0,2583 0,0176 0,0035 Quadratura Gaussiana Exemplo: Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 2 pontos. dxe x 10 0 16/05/2019 18 Quadratura Gaussiana Solução: a) Temos no intervalo [0,10]. Fazendo a mudança de variáveis xexf )( 1,1 temos10,0 para 55)( 2 1 )( 2 1 )( tx tabtabtx dtdtabdttxdx 5)( 2 1 )( )( 2 1 )( 2 1 )( 2 1 )( ababtabftF 5--5te5 Quadratura Gaussiana Solução: a) Para dois pontos: O erro verdadeiro: A Regra dos Trapézios necessitaria de 16 pontos para atingir este erro. Através da 1/3 de Simpson seriam necessários 9 pontos 3 3 3 3 )( 1 1 10 0 FFdttFIdxe Gauss x 5 3 3 55 3 3 5 55 ee 606102,0 393853,0606102,0999955,0Erro Quadratura Gaussiana Conclusões: As fórmulas da Quadratura Gaussiana produzem melhores resultados que aquelas dos métodos de Newton-Cotes com partição regulares (Trapézio, Simpson) Quando aumentamos o número de pontos todos métodos melhoram a precisão Se o intervalo for grande, no caso de Trapézio e Simpson Repetidas, podemos criar subintervalos e aplicar Quadratura Gaussiana em cada intervalo Problema: Se não tivermos f(x) e sim uma tabela de dados experimentais, então o método da Quadratura Gaussiana não é aplicável Quadratura Gaussiana Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 2 e 3 pontos. Resposta: 2 pontos 38,287 3 pontos 38,076 dxx 56 9 1 Quadratura Gaussiana Exercício: Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 2 e 3 pontos. Resposta: 2 pontos 1,718 3 pontos 1,7183 dxex 1 0 Quadratura Gaussiana - Exercícios 1) Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 2 pontos. Resposta: 0,746594 2) Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 2 pontos. Resposta: 27,14719 dxe x 21 0 dxxsenx x )( 4 2 0 16/05/2019 19 Quadratura Gaussiana – Exercícios 3) Considere a integral a) Calcule uma aproximação para a integral utilizando a Quadratura de Gauss para 3 pontos. Resposta: 3,14108 dx x2 1 0 1 1 4 Laboratorio_01_apres.pdf Cálculo Numérico Laboratório 1 Prof. Wellington Passos de Paula wpassos@gmail.com Introdução � Software destinado a fazer cálculos com matrizes ( MATLAB = MATrix LABoratory) � Comandos próximos à forma como escrevemos expressões � Torna o uso mais simples� Torna o uso mais simples Introdução � Tela Inicial do MATLAB � Desktop � Desktop Layout � Default Introdução � Janela de Comandos (Command Window) � Principal janela do MATLAB, onde os comandos devem ser executados � Utilizando as teclas e é possível navegar pelos comandos já digitados � Comando clc limpa a tela� Comando clc limpa a tela Introdução � Local de Trabalho (Workspace) � Armazena as variáveis já criadas, com seus respectivos valores � Comando clear (executado na Command Window) limpa todas as variáveis existentesvariáveis existentes � Histórico de Comandos (Command History) � Armazena os comandos já executados � Para limpar essa lista, botão direito � Clear Command History Introdução � Diretório de trabalho (Current Folder) � Se algum arquivo for chamado, este deve estar neste diretório � Atividade Inicial (Importante)� Atividade Inicial (Importante) � Fazer o apontamento do MATLAB para um diretório pessoal (preferencialmente no Desktop) � Sempre refazer esse apontamento no inicio de cada aula Cálculos Simples � MATLAB faz cálculos simples e científicos como uma calculadora � Ex: Suponha que você vá a uma loja e compre 3 objetos que custam 25 reais cada e 5 objetos que custam 12 reais cada. Quanto custou a sua compra? � MATLAB respeita a precedência de operadores e cria uma nova variável � Quando não definimos uma variável para guardar o resultado, o MATLAB utiliza a variável ans Cálculos Simples � Ex: Suponha que você vai a uma loja e compra 3 objetos que custam 25 reais cada e 5 objetos que custam 12 reais cada. Quanto custou a sua compra? Outra maneira de resolver o problema: Cálculos Simples � Utilizamos a virgula para separar os comandos na mesma linha � Poderíamos usar ponto e virgula também. Nesse caso, o MATLAB não mostra o resultado dos comandos � Quando definimos uma variável para receber o resultado de uma operação, o MATLAB armazena este resultado somente nessa variável, NÃO utilizando a variável ans Cálculos Simples � Para ver o valor contido em uma variável, verificamos o Workspace ou digitamos seu nome na linha de comando � Operações aritméticas disponibilizadas pelo MATLAB: + (soma), - (subtração), * (multiplicação), / (divisão), ^ (exponenciação) � Outras operações mais complexas, como cálculo de raiz, são disponibilizados como funções � Precedência de operações ^, * ou /, + ou – Nomenclatura de Variáveis � Nomes de variáveis devem ser iniciados por letras � Não podem conter espaços � Não podem conter caracteres de pontuação � MATLAB é case sensitive � casa ≠ cAsA � Variáveis pré-definidas: � ans, pi, eps, flops, inf, NaN, nan, nargin, nargout, realmin, realmax � MATLAB não impossibilita a redefinição destas variáveis, mas esta deva ser evitada Funções Científicas � Existem funções pré-definidas no MATLAB � Para se verificar a lista de todas as funções para manipulação de expressões algébricas, disponíveis no MATLAB, basta executar o comando help symbolic � Para se verificar a documentação de uma função específica, execute o comando help nome_funcao � Ex: help sqrt Funções Científicas � Alguns outros exemplos de funções científicas disponíveis: Formatos Numéricos � MATLAB respeita certas regras para apresentar resultados numéricos � Número Inteiro � Apresentado como Inteiro � Número Real � Apres. com aprox. de até quatro casas � Números Muito Grandes ou Muito Pequenos � Apresentados com Notação CientíficaApresentados com Notação Científica � Possível definição de novos formatos � format short – 5 dígitos (1 antes da virgula, + 4 decimais) � format long – 16 dígitos (1 antes da virgula, + 15 decimais) � format rat – formato racional � Formatos diferentes não alteram a representação interna dos números no MATLAB Programação no MATLAB � Não é preciso declarar o tipo das variáveis � O tipo é definido na atribuição do valor � x = 5 � x é uma variável inteira � y = 4.5 � y é uma variável do tipo doule � z = ‘Hello World’ � z é uma variável do tipo string � No MATLAB o primeiro índice de uma string é 1, ou seja z(1) retorna ‘H’ � Comando de atribuição = � Comando % indica comentário � MATLAB não aceita comandos como y--, x++, etc Programação no MATLAB � Estruturas de controle de fluxo � Desvio condicional Caso simples: Utilizando elseif: Com if aninhado: Exemplo: Programação no MATLAB � Estruturas de controle de fluxo � Repetição � For Caso simples: Exemplo: � No caso simples o incremento é de 1 em 1 Definindo o incremento: Exemplo: Programação no MATLAB � Estruturas de controle de fluxo � Repetição � While � O comando break pode ser utilizado para interromper o while durante sua execução (funciona para o for também) � Operadores relacionais: <, >, <=, >=, ==, ~= � Operadores lógicos: & (e), | (ou), ~ (não) Criação de Funções Externas � Executar muitos comandos na Command Window pode ser confuso � Solução: Criação de funções externas � File � New� Function Criação de Funções Externas � MATLAB já cria um esqueleto da função � output_args � lista de parâmetros de saida (pode ser um vetor ou uma matriz) � Untitled� nome sugerido para a função, pode ser mudado � Obs: Quando salvar o arquivo, utilizar o mesmo nome da � Obs: Quando salvar o arquivo, utilizar o mesmo nome da função criada) � input_args� lista de parâmetros de entrada � Os comentários nas linhas abaixo do titulo são mostrados quando o comando help nome_funcao é executado Criação de Funções Externas � Criando uma função que calcula a soma de x + y � Executando o help � Executando a função Laboratorio_02_apres.pdf Cálculo Numérico Laboratório 2 Prof. Wellington Passos de Paula wpassos@gmail.com Configurações Iniciais � Tela Inicial do MATLAB � Desktop � Desktop Layout � Default Equações não lineares - Gráficos � Comando para plotagem de gráficos: � plot(x, y) � comando para plotagem de gráficos � x e y podem ser dois vetores com os valores de cada ponto (x,y) no gráfico � Aplicação: problemas de interpolação (futuro) Equações não lineares - Gráficos � Comando para plotagem de gráficos: � plot(x, y) � comando para plotagem de gráficos � x pode ser um intervalo e y pode estar em função de x � x = -2:0.5:2 (executar sem o ; a fim de ver os resultados na tela) � Este comando gera uma sequência de números, iniciada em -2, com saltos de 0.5 em 0.5, até finalizar com o valor 2com saltos de 0.5 em 0.5, até finalizar com o valor 2 � y = x.^3-x-1 (executar sem o ; afim de ver os resultados na tela) � Este comando define o valor de y em função de x. Para cada valor de x gerado pelo primeiro comando, será gerado, respectivamente, um valor para y, a partir da formula acima. � Para os operadores * (multip), / (divisão) e ^ (exp.) , utilizar .*, ./ e .^ (operadores termo a termo) � manter essa prática sempre quando for plotar gráficos � grpx=plot(x,y); Equações não lineares - Gráficos � Comando para plotagem de gráficos: -2 0 2 4 6 y 3 4 � Uma vez que o gráfico esteja plotado, não fechar a janela � Necessário repetir os comandos -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 -8 -6 -4 -2 x1 2 3 4 50-1-2-3-4 1 2 3 -4 -3 -2 -1 Equações não lineares - Gráficos � Formatação de gráficos: � Justificativa para a necessidade de atribuir o resultado da função plot a uma variável (grpx=plot(x,y);) � Para alterar algumas propriedades, é necessário identificar o gráfico � Tamanho da linha, estilo e cor:� Tamanho da linha, estilo e cor: � set(grpx, 'LineWidth',3, 'LineStyle','--', 'Color','k'); � Alguns estilos de linha: '--' (tracejada) ':' (pontilhada) '-.' (tracejada-pontilhada) '-' (contínua) � Algumas cores 'k' (black) 'b' (blue) 'r' (red) 'g' (green) Equações não lineares - Gráficos � Formatação de gráficos: � Definições grid � axis([-4 5 -4 4]); � Limites dos eixos � grid ; � Mostrar as linhas do grid Equações não lineares - Gráficos � Formatação de gráficos: y 3 4 0 1 2 3 4 x1 2 3 4 50-1-2-3-4 1 2 3 -4 -3 -2 -1 -4 -3 -2 -1 0 1 2 3 4 5 -4 -3 -2 -1 Equações não lineares - Gráficos � Formatação de gráficos: � Legenda para a curva: � legend(grpx,'Curva teste'); � Título dos eixos: � xlabel('eixo X'); � ylabel('eixo Y'); � Título gráfico: � title('Gráfico teste'); Equações não lineares - Gráficos � Formatação de gráficos: y 3 4 0 1 2 3 4 e i x o y Gráfico teste Curva teste x1 2 3 4 50-1-2-3-4 1 2 3 -4 -3 -2 -1 -4 -3 -2 -1 0 1 2 3 4 5 -4 -3 -2 -1 eixo X e i x o y Equações não lineares - Gráficos � É possível fazer tudo de uma vez? � Para tipo e cor das linhas, sim: � plot(x, y, ':r'); 4 6 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 -8 -6 -4 -2 0 2 4 Equações não lineares - Gráficos � É possível fazer a edição de gráficos dentro da própria janela Equações não lineares - Gráficos � É possível colocar mais de uma curva dentro de um mesmo gráfico: � t = 0:0.1:4.*pi; z=sin(t); w=cos(t); g=0.*t; a = plot(t, z, t, w, t, g); legend(a, 'seno', 'cosseno', 'marcacao eixo X'); text(2.2, -0.5, 'cosseno'); 1 text(2.2, -0.5, 'cosseno'); 0 2 4 6 8 10 12 14 -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 cosseno seno cosseno marcacao eixo X Equações não lineares - Gráficos � É possível plotar gráficos em três dimensões � points = linspace(-2, 2, 40); [x, y] = meshgrid(points, points); z = 2./exp((x-.5).^2+y.^2)-2./exp((x+.5).^2+y.^2); surf(x, y, z); -2 -1 0 1 2 -2 -1 0 1 2 -1.5 -1 -0.5 0 0.5 1 1.5 Equações não lineares - Gráficos � Galeria de gráficos do MATLAB: � http://www.mathworks.com/discovery/gallery.html Equações não lineares - Raízes � O MATLAB possui duas poderosas funções para o cálculo de raízes de equações � Função roots ( p ) � Obtém todas as raízes de uma equação polinomial � Função fzero (f, x0) � Determina uma raiz de uma equação polinomial ou transcendente Equações não lineares - Raízes � Função roots (p) � Parâmetros � p ����O polinômio fornecido na entrada da função roots é representado através de um vetor de coeficientes, em suas respectivas posições � Dado p(x) = x4 + 2x3 - 13x2 - 14x + 24� Dado p(x) = x4 + 2x3 - 13x2 - 14x + 24 � p = [1 2 -13 -14 +24]; � r = roots(p) � Resultado mostrado pelo Matlab r = -4.0000 3.0000 -2.0000 1.0000 Equações não lineares - Raízes � Função roots (p) � Dado p(x) = x3 - x – 1 � p = [1 0 -1 -1]; � r = roots(p); � Resultado mostrado 3 4 Gráfico teste Curva teste pelo Matlab r = 1.3247 -0.6624 + 0.5623i -0.6624 - 0.5623i -4 -3 -2 -1 0 1 2 3 4 5 -4 -3 -2 -1 0 1 2 3 eixo X e i x o y Equações não lineares - Raízes � Função roots (p) � Dado p(x) = x3 - x – 1 � O comando poly(x) faz o caminho contrário, ou seja, recupera o polinômio a partir de suas raízes � p=poly(r) Resultado mostrado pelo Matlab� Resultado mostrado pelo Matlab p = 1.0000 0.0000 -1.0000 -1.0000 � A função roots é limitada, pois só funciona para funções polinomiais Equações não lineares - Raízes � Função fzero (f, x0) � Parâmetros � f���� Função de trabalho, fornecida como um string � Como a função é fornecida em um string, não é necessário, para os operadores ^ (exponenciação.) , * (multiplicação) e / (divisão) , utilizar o . (ponto)(divisão) , utilizar o . (ponto) � Mas e se eu quiser colocar o ponto, para padronizar? Sem problemas, não vai gerar erro. � x0���� Vetor aproximação da raiz. Pode conter 1 ou 2 elementos. Para o segundo caso (2 elementos), a raiz deve estar contida dentro do intervalo definido entre os dois elementos. Para o caso em que f tenha mais de uma raiz, os cálculos devem ser feitos separadamente Equações não lineares - Raízes � Função fzero (f, x0) � Dado f(x) = x3 - x – 1 � f='x^3-x-1'; � x=fzero(f, 1.5) � Resultado 3 4 Gráfico teste Curva teste mostrado pelo Matlab x = 1.3247 -4 -3 -2 -1 0 1 2 3 4 5 -4 -3 -2 -1 0 1 2 eixo X e i x o y Equações não lineares - Raízes � Função fzero � Dado f(x) = 3x2 * seno(x) * e-X � f='3*x^2*sin(x)*exp(-x)'; � x=fzero(f, 3.5) � Resultado 2 mostrado pelo Matlab x = 3.141592653589793 1 1.5 2 2.5 3 3.5 4 -1 -0.5 0 0.5 1 1.5 Equações não lineares - Raízes � Função fzero –Observações: � A função fzero funciona bem para funções continuas � Ex: fzero('tan', 1) � Resultado mostrado pelo Matlab ans = 1.570796326794898 � A função tangente não é continua em π/2 � A função fzero necessita de troca de sinal da função para encontrar a raiz � Ex: f(x) = x2 = 0 � x = 0 � y =fzero('x^2', 0.5) � Resultado mostrado pelo Matlab y = NaN Equações não lineares - Raízes � Função fzero (f, x0) - Passos na execução da função � Passo 1: Análise crítica � Parece existir algum valor de x tal que f(x) = 0? � Sim/Talvez � Passo 2 � Não � A função não tem raiz � Passo 2: Plotar o gráfico � Importante: � Ao definir o intervalo, verificar se a função é contínua, para aquele intervalo � Verificar se a curva cruza o eixo x, pois caso contrário, a função fzero (f, x0) falha � A partir do gráfico plotado, definir um valor inicial (chute) x0 � Passo 3: Executar a função fzero (f, x0) Laboratorio_02_FundamentosParaAnaliseFuncoes.pdf 06/09/2018 1 Fundamentos para análise de funções no Octave/Matlab Funções e gráficos • O gráfico de uma função é gerado à partir dos seus pontos (xi, yi); – Define-se o conjunto de pontos xi: • x = -2:0.5:2; – Este comando gera uma sequência de números, iniciada em -2, com saltos de 0.5 em 0.5, até finalizar com o valor 2. – Define-se o conjunto de pontos yi correspondentes, utilizando a função em questão. • Para a função f(x)=x3-x-1: – y = x.^3-x-1; » Este comando define o valor de y em função de x. Para cada valor de x gerado pelo primeiro comando, será gerado, respectivamente, um valor para y, a partir da formula acima. » Para os operadores * (multip), / (divisão) e ^ (exp.) , utilizar .*, ./ e .^ (operadores termo a termo) Funções e gráficos • Os vetores x e y contem todos os pontos do gráfico; • Para plotar, use o comando plot(x,y) Funções e gráficos • Para manipular as propriedades do gráfico, gere-o assim: grafico = plot(x,y); • A variável grafico permitirá realizar diversas alterações: set(grafico, 'LineWidth',3, 'LineStyle','--', 'Color','k'); • Alguns estilos de linha: '--' (tracejada) ':' (pontilhada)'-.' (tracejada-pontilhada) '-' (contínua) • Algumas cores 'k' (black) 'b' (blue) 'r' (red) 'g' (green) Funções e gráficos axis([-4 5 -4 4]); % Limites dos eixos grid;% Mostrar as linhas do grid Funções e gráficos • Legenda para a curva: legend(grafico,'Curva teste'); • Título dos eixos: xlabel('eixo X'); ylabel('eixo Y'); • Título gráfico: title('Gráfico teste'); 06/09/2018 2 Funções e gráficos • É possível colocar mais de uma curva dentro de um mesmo gráfico: t = 0:0.1:4.*pi; z=sin(t); w=cos(t); g=0.*t; a = plot(t, z, t, w, t, g); legend(a, 'seno', 'cosseno', 'marcacao eixo X'); text(2.2, -0.5, 'cosseno'); grid; Exercícios práticos • Plote as seguintes funções: • Dica: é preciso definir o intervalo numérico em x, gerando um vetor, para depois lançar estes números num vetor y, correspondendo a cada uma das funções. Manipulação de funções • Podemos definir uma função como simbólica (para integração, derivação, etc.): – Depende do pacote “symbolic”. syms f(x); % define f(x) como simbólica f(x) = sin(x^2); % define a função propriamente dita df = diff(f,x); % calcula a derivada em relação à variável x df(1) % a derivada recebe o valor 1 como parâmetro Manipulação de funções • Criando uma função f(x): function retval = polinomio(x) retval = x.^2+x-6; endfunction • Informando ao Octave que polinomio é uma função manipulável, criando o ponteiro f: f=@polinomio Manipulação de funções • Como f faz referência à função polinomio, podemos utilizá-la em outra função: function retval = bissecao (a, b, f, precisao) if(f(a)*f(b)>0) fprintf('Intervalo invalido (f(a)*f(b)>0)! %f * %f = %f > 0\n', f(a), f(b), f(a)*f(b)); return; end x = (b+a)/2; k = 0; fprintf('Iteracao | \t x \t\t | \t f(x) \t\t | \t |b-a| \t\t \n'); while( abs(b-a) >= precisao && abs(f(x)) >= precisao ) if(f(a)*f(x) < 0) % trocar o b b = x; else % trocar o a a = x; end x = (b+a)/2; k = k+ 1; fprintf(' %d \t | %.5f \t\t | %.5f \t\t | %.5f \t\t \n', k, x, f(x), abs(b-a)); end retval = x; fprintf('resposta x = %f, pois f(%f) = %f\n', x, x, f(x)); endfunction Laboratorio_03_Obtencao de raizes.pdf 06/09/2018 1 Obtenção de raízes Raízes usando Octave • Conforme visto, a análise do gráfico da função f(x) e demais (ф(x), ф’(x), f’(x), f’’(x)) é relevante para obtenção da raiz da função. – Gera-se um vetor de valores x e depois os valores y correspondentes. • Podemos implementar nossas próprias funções para obtenção de raízes. • Ou reusar funções existentes. Implementando um método... • Pode-se implementar um método para obtenção de raízes e aplicá- lo na função desejada. • Gera-se um arquivo contendo a função: function retval = polinomio(x) retval = x.^2+x-6; endfunction • Cria-se um ponteiro para a função: f=@polinomio • Implementa-se o método (slides a seguir). Bisseção function retval = bissecao (a, b, f, precisao) if(f(a)*f(b)>0) fprintf('Intervalo invalido (f(a)*f(b)>0)! %f * %f = %f > 0\n', f(a), f(b), f(a)*f(b)); return; end x = (b+a)/2; k = 0; fprintf('Iteracao | \t x \t\t | \t f(x) \t\t | \t |b-a| \t\t \n'); while( abs(b-a) >= precisao && abs(f(x)) >= precisao ) if(f(a)*f(x) < 0) % trocar o b b = x; else % trocar o a a = x; end x = (b+a)/2; k = k+ 1; fprintf(' %d \t | %.5f \t\t | %.5f \t\t | %.5f \t\t \n', k, x, f(x), abs(b-a)); end retval = x; fprintf('resposta x = %f, pois f(%f) = %f\n', x, x, f(x)); endfunction Posição Falsa function retval = posicaoFalsa (a, b, f, precisao) if(f(a)*f(b)>0) fprintf('Intervalo invalido (f(a)*f(b)>0)! %f * %f = %f > 0\n', f(a), f(b), f(a)*f(b)); return; end x = (a*f(b)-b*f(a))/(f(b)-f(a)); k = 0; fprintf('Iteracao | \t x \t\t | \t f(x) \t\t | \t |b-a| \t\t \n'); while( abs(b-a) >= precisao && abs(f(x)) >= precisao ) if(f(a)*f(x) < 0) % trocar o b b = x; else % trocar o a a = x; end x = (a*f(b)-b*f(a))/(f(b)-f(a)); k = k+ 1; fprintf(' %d \t | %.5f \t\t | %.5f \t\t | %.5f \t\t \n', k, x, f(x), abs(b-a)); end retval = x; fprintf('resposta x = %f, pois f(%f) = %f\n', x, x, f(x)); endfunction Newton-Raphson function retval = NewtonRaphson (x, f, precisao) pkg load optim; % para uso da função deriv x = x-f(x)/deriv(f, x); k = 0; fprintf('Iteracao | \t x \t\t | \t f(x) \t\t \n'); fprintf(' %d \t | %.5f \t\t | %.5f \t\t \n', k, x, f(x)); while( abs(f(x)) >= precisao ) x = x-f(x)/deriv(f, x); k = k+ 1; fprintf(' %d \t | %.5f \t\t | %.5f \t\t \n', k, x, f(x)); end retval = x; fprintf('resposta x = %f, pois f(%f) = %f\n', x, x, f(x)); endfunction 06/09/2018 2 Secante function retval = secante (x_k_1, x_k, f, precisao) x = (x_k_1*f(x_k)-x_k*f(x_k_1))/(f(x_k)-f(x_k_1)) k = 0; fprintf('Iteracao | \t x \t\t | \t f(x) \t\t \n'); while(abs(x_k - x_k_1) >= precisao && abs(f(x)) >= precisao ) x = (x_k_1*f(x_k)-x_k*f(x_k_1))/(f(x_k)-f(x_k_1)); x_k_1 = x_k; x_k = x; k = k+ 1; fprintf(' %d \t | %.5f \t\t | %.5f \t\t \n', k, x, f(x)); end retval = x; fprintf('resposta x = %f, pois f(%f) = %f\n', x, x, f(x)); endfunction Funções oferecidas pelo Octave/Matlab • Função roots ( p ) – Obtém todas as raízes de uma equação polinomial • Função fzero (f, x0) – Determina uma raiz de uma equação polinomial ou transcendente Funções oferecidas pelo Octave/Matlab • Função roots (p) – Parâmetros • p -> O polinômio fornecido na entrada da função roots é representado através de um vetor de coeficientes, em suas respectivas posições – Dado p(x) = x4 + 2x3 - 13x2 - 14x + 24 p = [1 2 -13 -14 +24]; r = roots(p) • Resultado mostrado pelo Matlab r = -4.0000 3.0000 -2.0000 1.0000 Funções oferecidas pelo Octave/Matlab • Função roots (p) – Dado p(x) = x3-x-1 p =[1 0 -1 -1]; r = roots(p) r = 1.32472 + 0.00000i -0.66236 + 0.56228i -0.66236 - 0.56228i Funções oferecidas pelo Octave/Matlab Funções oferecidas pelo Octave/Matlab • Função fzero(f,x0) – Parâmetros: • f -> Função de trabalho fornecida como uma handle function: f =@(x) x.^3-x-1; f =@(x) 2*sin(x); • x0 -> Ponto inicial: >> f = @(x) x.^2 f = @(x) x .^ 2 >> fzero(f, 1) ans = 0 >> 06/09/2018 3 Funções oferecidas pelo Octave/Matlab • Função fzero(f,x0) – Dado f(x)=x3-x-1 f =@(x) x.^3-x-1; x=fzero(f,1.5) x = 1.3247 Funções oferecidas pelo Octave/Matlab • Função fzero(f,x0) – Dado f(x)=3x2*seno(x)*e-x f =@(x) 3*x.^2*sin(x)*exp(-x); x=fzero(f,3.5); x = 3.14159 Funções oferecidas pelo Octave/Matlab Laboratorio_04_apres.pdf Cálculo Numérico Laboratório 3 Prof. Wellington Passos de Paula wpassos@gmail.com Configurações Iniciais � Tela Inicial do MATLAB � Desktop � Desktop Layout � Default Configurações Iniciais � Diretório de trabalho (Current Folder) � Se algum arquivo for chamado, este deve estar neste diretório � Atividade Inicial (Importante) � Fazer o apontamento do MATLAB para um diretório pessoal (preferencialmente no Desktop) � Sempre refazer esse apontamento no inicio de cada aula Matrizes no Matlab � Matrizes são tratadas da mesma forma que escalares e vetores no Matlab, já que os últimos são casos particulares de matrizes � Criação de matrizes no Matlab � . � Ex: a=[1 2 3; 4 5 6] � É também possível criar matrizes digitando cada linha separadamente (utilizar Shift+Enter para mudar de linha) � Ex: b=[7 8 9 10 11 12] ];;[ 212222111211 mnmmnn aaaaaaaaaNomeMatriz LLL= Matrizes no Matlab � Outros comandos para criação de matrizes: � eye(n) � cria uma matriz identidade de ordem n � ones(n)� cria uma matriz de ordem n com todos os elementos iguais a 1 � zeros(n)� cria uma matriz nula de ordem n � rand(n)� matriz aleatória (entre 0 e 1) de ordem n � magic(n)� cria uma matriz aleatória ordem n. O somatório de cada linha é o mesmo de cada coluna Matrizes no Matlab � O comando size permite obter as dimensões da matriz � Ex: [l,c]=size(b) � Retorno do Matlab: l = 2 c = 3 � Para acessar uma determinada posição da matriz, basta digitar o comando NomeMatriz(i,j), lembrando que, no Matlab, os índices começam a partir de 1 � Ex: b(2,1) � Retorno do Matlab: ans = 10 Matrizes no Matlab � Para retornar de uma vez todos os elementos da linha/coluna da matriz: � Ex: b(1,:) ���� retorna os elementos da primeira linha de b � Retorno do Matlab: ans = 7 8 9 � Ex: b(:,2) ���� retorna os elementos da segunda coluna de b � Retorno do Matlab: ans = 7 10 � É possível também remover linhas/colunas da matriz � Ex: c = b(:,[1 3]) � remove a segunda coluna da matriz b, atribuindo o resultado à variável c. � Retorno do Matlab: c = 7 9 10 12 Matrizes no Matlab - Operações � As operações matriciais são executadas de maneira semelhante às operações escalares � Ex: a=[1 2; 3 4]; c=[3 5; -5 2]; � Adição: d = a + c � Retorno do Matlab: d = 4 7 -2 6 � Subtração: d = c - a � Retorno do Matlab: d = 2 3 -8 -2 Matrizes no Matlab - Operações � As operações matriciais são executadas de maneira semelhante às operações escalares � Ex: a=[1 2; 3 4]; c=[3 5; -5 2]; � Multiplicação: d = a * c � Retorno do Matlab: d = -7 9 -11 23 E se fizermos d = a .* c, o que esperar? � Divisão: d = c/a � Retorno do Matlab: d = 1.5000 0.5000 13.0000 -6.0000 Matrizes no Matlab - Operações � As operações matriciais são executadas de maneira semelhante às operações escalares � Ex: a=[1 2; 3 4]; c=[3 5; -5 2]; � Exponenciação: d = a ^ 2 � Retorno do Matlab: d = 7 10 15 22 � E se fizermos a .^ 2, o que esperar? � Transposta: d=d' � Retorno do Matlab: d = 7 15 10 22 Matrizes no Matlab - Operações � As operações matriciais são executadas de maneira semelhante às operações escalares � Ex: a=[1 2; 3 4]; c=[3 5; -5 2]; � Inversa: d = inv(a) � Retorno do Matlab: d = -2.0000 1.0000 1.5000 -0.5000 � Determinante: x=det(a) � Retorno do Matlab: x = -2 Solução de Sistemas Lineares � O Matlab dispõe de operadores e funções que podem ser utilizadas para solucionar sistemas de equações lineares � Divisão à esquerda � Permite resolver um sistema linear de n equações, no formato Ax = b, onde: A = matriz de coeficientes n x n x = matriz coluna (n x 1) de incógnitas b = matriz coluna (n x 1) de termos independentes � O comando da divisão à esquerda tem a seguinte forma: x = a\b � A solução é calculada fazendo-se a Eliminação de Gauss Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações: = × 6 5 43 21 2 1 x x Solução: a= [1 2; 3 4]; b=[5 ; 6] x=a\b Retorno do Matlab: x = -4.0000 4.5000 Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações: − = × −− −− −− 17 16 5.6 12 15.152212 5.525.65.71 65.676 6324 4 3 2 1 x x x x Solução: a=[4 -2 -3 6; -6 7 6.5 -6; 1 7.5 6.25 5.5;-12 22 15.5 -1] b=[12 ; -6.5; 16; 17] x=a\b Retorno do Matlab: x = 2.0000 4.0000 -3.0000 0.5000 −− 1715.152212 4x Solução de Sistemas Lineares � Decomposição LU com Pivoteamento Parcial � Comando para decompor a matriz A em L e U [L,U,P] = lu(A) onde: L = matriz triangular inferior unitária U = matriz triangular superior P = matriz identidade reordenada (pivoteamento parcial) � A partir do cálculo acima, a resolução dos sistemas Ly=Pb e Ux=y é realizada através da divisão à esquerda. Solução de Sistemas Lineares = × 6 5 43 21 2 1 x x � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: Solução: a= [1 2; 3 4]; b=[5 ; 6]; [l,u,p] = lu(a) Retorno do Matlab: l = 1.0000 0 u= 3.0000 4.0000 p= 0 1 0.3333 1.0000 0 0.6667 1 0 Solução de Sistemas Lineares = × 6 5 43 21 2 1 x x � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: Solução: Resolvendo o sistema Ly=Pb y= l\(p*b) Retorno do Matlab: y = 6 3 Solução de Sistemas Lineares = × 6 5 43 21 2 1 x x � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: Solução: Resolvendo o sistema Ux=y x= u\y Retorno do Matlab: x= -4.0000 -4.5000 Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: − = × −− −− 16 5.6 12 5.525.65.71 65.676 6324 3 2 1 x x x Solução: a=[4 -2 -3 6; -6 7 6.5 -6; 1 7.5 6.25 5.5;-12 22 15.5 -1] b=[12 ; -6.5; 16; 17] −− 1715.152212 4 3 x Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: − = × −− −− 16 5.6 12 5.525.65.71 65.676 6324 3 2 1 x x x Solução: [l,u,p] = lu(a) Retorno do Matlab: −− 1715.152212 4 3 x Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: − = × −− −− 16 5.6 12 5.525.65.71 65.676 6324 3 2 1 x x x Solução: Resolvendo o sistema Ly=Pb y= l\(p*b) Retorno do Matlab: y = 17.0000 17.4167 7.7143 -0.4000 −− 1715.152212 4 3 x Solução de Sistemas Lineares � Exemplo: Encontre as raízes do sistema de equações utilizando a decomposição LU: − = × −− −− 16 5.6 12 5.525.65.71 65.676 6324 3 2 1 x x x Solução: Resolvendo o sistema Ux=y x= u\y Retorno do Matlab: x = 2.0000 4.0000 -3.0000 0.5000 −− 1715.152212 4 3 x Laboratorio_05_apres.pdf Cálculo Numérico Laboratório 4 Prof. Wellington Passos de Paula wpassos@gmail.com Configurações Iniciais � Tela Inicial do MATLAB � Desktop � Desktop Layout � Default Interpolação � O MATLAB disponibiliza métodos para interpolação de funções � Neste laboratório, trabalharemos com as funções abaixo: � polyfit (x, y, m) � Calcula o polinômio de grau até n-1 que passa pelo conjunto de n pontos (x, y) � interp1(x, y, xi, 'método') � Dado um conjunto de pontos (x,y),faz a interpolação dos pontos contidos na variável xi, utilizando para isso o método de interpolação informado em 'método' Interpolação - polyfit � Função polyfit (x, y, m) � Parâmetros � x� vetor contendo as coordenadas horizontais dos pontos � y� vetor contendo as coordenadas verticais dos pontos � m� grau do polinômio interpolado Interpolação - polyfit � Função polyfit (x, y, m) � Observação: � Embora o polinômio interpolado possa ter grau até n-1, onde n é o número de pontos fornecidos, normalmente utilizar polinômios de grau muito elevados não é uma boa escolha ( problemas de arredondamento, etc ). Assim, escolha ( problemas de arredondamento, etc ). Assim, costumamos trabalhar com polinômios de graus mais baixos, como 3,4 e 5 que, na maioria dos casos, apresentam boas aproximações para a função interpolada. Interpolação - polyfit � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: x 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 y 0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 Solução: Gerando os valores de x e y: x = 0:0.4:3.6 y = [0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6] Gerando o polinômio interpolador: p =polyfit(x, y, 4) Interpolação - polyfit � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: x 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 y 0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 Solução: Resultado do MATLAB: p = 0.2726 -0.5966 -2.4014 7.9459 0.0591 O resultado acima mostra que, para os pontos (x,y) fornecidos, o MATLAB calculou o seguinte polinômio de grau 4: P4(x) = 0.2726x4 - 0.5966x3 - 2.4014x2 + 7.9459x + 0.0591 Interpolação - polyfit � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: x 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 y 0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 Solução: Agora, vamos analisar a qualidade do polinômio retornado pelo MATLAB. Para isso, vamos plotar, num mesmo gráfico, o pontos (x,y) fornecidos na entrada e os pontos (x, yi), onde yi é calculado a partir da avaliação do polinômio P4 no ponto x. Interpolação - polyfit � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: x 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 y 0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 Solução: Vamos utilizar a função polyval (p, x), onde p é o polinômio que vamos avaliar e x é o conjunto de 1 ou mais pontos nos quais desejamos avaliar p. Logo: yi = polyval (p, x) Resultado do Matlab: yi = 0.0591 2.8220 4.6851 5.6706 5.9679 5.9344 6.0945 7.1403 9.9314 15.4948 Interpolação - polyfit � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: x 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 y 0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 Solução: Plotando o gráfico: plot(x, y, '*', x, yi, '-') No gráfico, os pontos (x, y) originais estarão marcados com *, já os pontos (x, yi), calculados pelo polinômio interpolado, estarão ligados por uma linha contínua Interpolação - polyfit 10 12 14 16 � Considere o conjunto de pontos abaixo e obtenha o polinômio interpolador de grau 4: 0 0.5 1 1.5 2 2.5 3 3.5 4 0 2 4 6 8 10 E aí, o polinômio interpolado está bom ou ruim? Interpolação – interp1 � Função interp1(x, y, xi, 'método') � Parâmetros � x� vetor contendo as coordenadas horizontais dos pontos � y� vetor contendo as coordenadas verticais dos pontos � xi� conjunto de 1 ou mais pontos cujos quais desejamos interpolar � 'método' � método utilizado na interpolação. Existem vários métodos disponíveis, mas nós vamos trabalhar com somente 2: � 'linear' � interpolação spline linear � 'cubic' � interpolação spline cúbica Interpolação – interp1 � Função interp1(x, y, xi, 'método') � Observações: � O vetor x deve ser monotônico ( os elementos devem estar na ordem ascendente ou descendente ) Interpolação – interp1 � Considere o conjunto de pontos abaixo e realize a interpolação para o ponto x = 18 x 8 11 14 17 20 y 5 9 10 8 7 Solução: Gerando os valores de x e y: x = 8:3:20 y = [5 9 10 8 7] Interpolação – interp1 � Considere o conjunto de pontos abaixo e realize a interpolação para o ponto x = 18 x 8 11 14 17 20 y 5 9 10 8 7 Solução: Calculando a interpolação no ponto x = 18 yi = interp1(x,y,18,'linear') ���� método linear Resposta do MATLAB: yi = 7.6667 yi = interp1(x,y,18,'cubic') ���� método cúbico Resposta do MATLAB: yi = 7.5802