Baixe o app para aproveitar ainda mais
Prévia do material em texto
PONTIFÍCIA UNIVERSIDADE CATÓLICA DO RIO GRANDE DO SUL FACULDADE DE MATEMÁTICA CÁLCULO NUMÉRICO Notas de Aula – Aplicações – Exercícios Eliete Biasotto Hauser Índice 1 Sistema de Ponto Flutuante Normalizado – Teoria dos Erros....................4 2 Resolução de Equações Algébricas e Transcendentes.............................. 9 3 Sistemas de Equações Lineares................................................................ 25 4 Interpolação Polinomial............................................................................ 36 5 Ajuste de Funções..................................................................................... 46 6 Integração Numérica.................................................................................60 5 Resolução Numérica de Equações Diferenciais Ordinárias..................... 65 8 Resolução Numérica de Equações Diferenciais Parciais.......................... 72 9 Bibliografia ..............................................................................................85 Formulário................................................................................................... 86 Laboratórios utilizando Sistema Maple....................................................... 90 1 –Teoria dos Erros - Sistema de Ponto Flutuante Um número y na base 2≥β , y = LL 2b1b,0a1a2a3a1nana − pode ser descrito na forma y = LL +−+−++++++−−+ 22b11b0a1a22a33a1n1nanna βββββββ . Por exemplo, 2106110271012105310326,3517 −×+−×++×+×+×= Em aritmética de ponto flutuante normalizado de t dígitos, y tem a forma: y = etd2d1d,0 ββ ×⎟⎠ ⎞⎜⎝ ⎛ L i) td2d1d,0 L é a mantissa (uma fração na base β), ii) 01d,1jd0 ≠−≤≤ β , j=1,2,...,t iii) “e” é um expoente inteiro que varia no intervalo [m,M]. M e m dependem da máquina utilizada. Em geral, m = -M∈Z. iv) t define a precisão da máquina, é o número de dígitos da mantissa. Obs: precisão ≠ exatidão(depende da precisão da máquina e do algorítmo utilizado) A união de todos os números de ponto flutuante normalizados com o zero: m vezest 0000.00 β×= 43421 L é chamado Sistema de Ponto Flutuante Normalizado e representado por F(β, t, m, M). Alguns exemplos de máquinas com precisão simples: a) HP 48 : F(10, 12, -498, 500) b) IBM 3090 : F(16, 6, -64, 63) c) Cray1 : F(2, 48, -8192, 8191) d) Burroughs B6700: F(8, 13, -63, 64) Em F valem as propriedades: 1) 0,1 x mβ é o menor número em módulo, não nulo, de F. 2) 0, 444 3444 21 vezest )1)...(1)(1( −−− βββ x Mβ é o maior número, em módulo, não nulo, de F. 3) # F = 1)1(1)1(2 +⎥⎦ ⎤⎢⎣ ⎡ +−−− mMtββ é o cardinal de F 4) Para qualquer mantissa, vale 1mantissa1 <≤−β . 5) Se Fy∈ , então Fy∈− . 6) F1eF0 ∈∈ . Se o expoente da base não pertencer a [m,M], y não pode ser representado em F. São os casos de erro de: - underflow, se e < m (ultrapassa a capacidade mínima) - overflow, se e > M (ultrapassa a capacidade máxima) 4 E.B.Hauser – Cálculo Numérico Se a representação do real y em F não é exata, é necessário utilizar um arredondamento. Os tipos de arredondamento mais conhecidos são: - Arredondamento para número mais próximo de máquina (Oy). - Arredondamento por falta, ou truncamento (∇ y). - Arredondamento por excesso (∆ y) Em F, geralmente, as operações de adição e multiplicação não são comutativas, associativas e nem distributivas, pois numa série de operações aritméticas, o arredondamento é feito após cada operação. Ou seja, nem sempre as operações aritméticas válidas para os números reais são válidas em F. Isto influi na solução obtida através de um método numérico. Assim, métodos numéricos matematicamente equivalentes, podem fornecer resultados diferentes. Medidas de Exatidão Quando se aproxima um número real x por x*, o erro que resulta é x-x*. Define-se: erro absoluto: EA = *xx − e o erro relativo: ER = x xx *− para 0x ≠ . A fim de ver o tipo de situação que pode ocorrer um erro relativo de grande magnitude, vamos considerar a diferença entre os números a seguir, por exemplo: x = 0,3721478693 y = 0,3720230572 yx − = 0,0001248121=0,1248121 x 10-3 Se os cálculos forem feitos em F(10, 5, -499, 499) com arredondamento Ox: x*= 0,37215, y*= 0,37202 e x*-y*= 0,00013 = 0,13000 x 10-3 Assim, o erro relativo entre os dois resultados é grande: %4 3-10 x 0,1248121 3-10 x 13000,03-10 x 0,1248121 yx -y*)*x(yx ≈−=− −− Na resolução de um problema o valor exato da solução x pode ser desconhecido. Podemos usar duas aproximações sucessivas de x, definindo: ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ + −+++−=+ 1 1log3,0)1,( kx kxkx kxkxDIGSE µ o qual expressa o número de dígitos significativos exatos de kx em relação a 1kx + . Aqui µ representa a unidade de arredondamento da máquina (µ = t1 2 1 −β se o arredondamen- to for Ox ). 5 1 - Teoria dos Erros Algoritmos Numéricos O Cálculo Numérico tem por objetivo estudar e aplicar algoritmos numéricos para a solução de problemas, visando o menor "custo" e confiabilidade do resultado (observar tempo de execução, número de operações aritméticas, quantidade de memória, propagação do erro, etc). Um bom algoritmo numérico deve se caracterizado por: a) Independência da máquina :a execução do programa pode ser realizada em diferentes máquinas. b) Inexistência de Erro Aritmético: problemas de overflow e underflow devem ser detecta- dos a priori c) Inexistência de Erro Lógico: previsão de tudo o que poderá ocorrer durante o processo (divisão por zero, por exemplo) d) Quantidade finita de cálculos. e) Existência de um critério de exatidão. f) O erro deve convergir para zero com precisão infinita. g) Eficiência: produz a resposta correta com economia Passos para a resolução de um problema: tipos de erros A resolução de problemas reais envolve diversas etapas: Estudaremos algoritmos numéricos a fim de obter uma solução numérica do problema, a qual, geralmente, difere da solução exata. Possíveis fontes de erro que geram essa diferença são: a) Simplificação do modelo matemático ( algumas variáveis envolvidas são desprezadas); b) Erro nos dados de entrada ( com frequência provindos de dados experimentais); c) Truncamento (por exemplo, substituição de um processo infinito por um finito e linearizações); d) Erro de arredondamento em aritmética de ponto flutuante. Curiosidades: Some disasters caused by numerical errors. http://ta.twi.tudelft.nl/users/vuik/wi211/disasters.html Aplicação O volume de uma esfera de raio r pode ser obtido de 3 3r4 V π= . Estimar o volume do hemisfério de raio “e” representado graficamente ao lado. Utilizar arredondamento para número mais próximo de máquina (0x) em F( 10 , 4 , -98 , 98). Problema Real Modelagem Matemática Solução 6 E.B.Hauser – Cálculo Numérico Exercícios 1) No sistema MapleV estimar 3.8e− utilizando ∑∞ = −=− 0i !i i)x(xe e ∑∞ = ==− 0i !i ix 1 xe 1xe com 26 termos cada e comparar com 3.8e− ≈0.2485168271x10-3. 2) Em F(10 , 3 , -98 , 98) e arredondamento por truncamento estimar p(2.73) se: a) 55.0x62x53x)x(p ++−= b) 55.0x)6x)5x(()x(p ++−= Em ambos os casos estimar o erro absoluto ao comparar com p(2.73) ≈0.11917x10-1. Obs: Estimar p(x) pelo algoritmo usual 01 2 2 3 3... 1 1)( axaxaxa nxna nxnaxp +++++−−+= exige n adiçõese (n2+n)/2 multiplicações enquanto que o algoritmo de Horner { 0)1)2...)2)1 1 (((...()( axaxaxnaxnaxna n xp +++−+−+ − = requer n adições e n multiplicações. 3) Sejam A, B, C e D matrizes genéricas de ordem 10x20, 20x50, 50x1 e 1x100 respectivamente. Utilizando a propriedade associativa, pode-se determinar o produto matricial AxBxCxD de diversas formas. Qual das duas abaixo é mais eficiente? Porque? a) Ax(Bx(CxD)) b) (Ax(BxC))xD 4) Representar o número real x na base 2 usando 8 algarismos significativos? Essa representação é exata? a) x=0.6 b) x=13.25 c) x= 2.47 5) Determinar o cardinal , regiões de underflow e overflow e todos elementos reais de: a) F(2,3,-1,2) b) F(3,2,-1,2) c) F(2,2,-2,2) 6) Representar, se possível, os números abaixo em utilizando arredondamento por truncamento( x∇ )e arredondamento para número mais próximo de máquina (Ox) em F(10,5,-2,2). a) 3 b) 3 200 c) 3 2000 d) e e) 3000 1 f) 2 7) Considerando: ∑ = = 10 12 1 i i A a) Calcular o valor de A utilizando precisão infinita.. b) Utilizando arredondamento por truncamento ( x∇ ) em F(10,3,-98,98), estimar o valor de A somando da direita para esquerda e após somando da esquerda para a direita. Comparar os resultados. 7 1 - Teoria dos Erros Respostas: 1) exp(-8.3); .0002485168271 > f1:=sum(((-x)^i)/i!, i=0..25): > f1a:=unapply(f1,x): > f1a(8.3); -.001241405905 Obs: Causas desse erro: subtração de grandezas muito próximas e adição de grandezas de diferentes ordens. > > f2:=1/(sum(((x)^i)/i!, i=0..25)): > f2a:=unapply(f2,x): > f2a(8.3); .0002485170000 2) a)p(2.73)= -0.05 ,erro absoluto = 0.061917 b)p(2.73)=0.032 ,erro absoluto = 0.020083 3) (Ax(BxC))xD é mais eficiente pois exige 2200 multiplicações enquanto que para calcular o produto Ax(Bx(CxD)) são necessárias 125000 multiplicações . OBS: Se M é de ordem pxq e N de ordem qxr, então MxN, de ordem pxr, é obtida efetuando pqr operações de multiplicações de elementos de M e N. 4-a) ( , ) ( , )0 6 10 0 10011001 2≈ b) 210 )01,1101()25,13( ≈ c) ( . ) ( , )2 47 10 10 011110 2≈ 5)-a) 0, 1/4, 1/2, 1, 2, 5/16, 5/8, 5/4, 5/2, 3/8, 3/4, 3/2, 3, 7/16, 7/8, 7/4, 7/2 e simétricos. # F = 33 Região de oferflow: ( , / ) ( / , )−∞ − ∪ +∞7 2 7 2 Região de underflow: (-1/4,1/4) - {0} b) 0, 1/9, 1/3, 1,3, 4/27, 4/9,4/3, 4, 5/27, 5/9, 5/3, 5, 2/9, 2/3, 2, 6, 7/27, 7/9, 7/3, 7, 8/27, 8/9, 8/3, 8 e seus simétricos. # F = 49 Região de oferflow: ( , ) (8, )−∞ − ∪ +∞8 Região de underflow: (-1/9,1/9) - {0} c) 0, 1/8, 1/4, 1/2, 1, 2, 3/16, 5/16, 3/4, 3/2, 3, e seus simétricos. # F = 21 Região de oferflow: ( , ) ( , )−∞ ∪ +∞3 3 Região de underflow: (-1/8,1/8) - {0} 6-a) 0,17320*101 e 0,17321 *101 b) 0,66666*102 e 0,66667*102 c) overflow (0,66666*103 ∉ e 0,66667*103 ∉ F) d) 0.27182*101 e 0.27183*102 e) underflow (0.33333*10 –3 ∉ F) f) 0,14142*101 e 0,14142 *101 7-a) 0,9990234375 b) 0,999 e 0,998 8 2. Resolução de Equações Algébricas e Transcendentes Objetivo: Resolver equações de forma ( ) 0=xf , isto é, determinar os valores de x para os quais a igualdade ( ) 0=xf é satisfeita. Se a função ( )xf só contém operações algébricas repetidas um número finito de vezes, a equação é dita algébrica. Ex.: polinômios, 0x x 1 ,01xx 53 2 =−=−− Equações Transcendentes são aquelas em que a variável independente x está submetida a uma operação não algébrica um número finito de vezes. Ex.: 0esenx ,0tgxxln ,01xe x2x =−=+=++ 2.1- Equações Polinomiais Revisão sobre polinômio: Seja ( ) 012211 axaxaxaxaxp nnnn +++++= −− K um polinômio de grau n, onde os coeficientes , são números reais e 0≠na . Temos que: 1) ( )xp possui n raízes, contadas as multiplicidades. 2) Se nxxx ,,, 21 K forem raízes reais de ( )xp , então ( )xp pode ser fatorado na forma: ( ) ( )( ) ( )n21n xxxxxxaxp −−−= K . Ex1: ( ) 2xx2xxp 23 −−+= , tem raízes: 2 ,1 ,1 321 −=−== xxx . Logo, ( ) ( )( )( )2x1x1xxp ++−= 3) Se bia + é raiz de ( )xp então biaz −= também o é. Ex.: ( ) 10x6xxp 2 +−= tem raízes i±3 . 9 E.B.Hauser – Cálculo Numérico 4) Se bia + é raiz de ( )xp de grau 2≥n , então ( )xp pode ser fatorado: ( ) ( ) ( )xqbaax2xxp 222 ⋅++−= onde o grau de ( )xq é 2−n . Ex.: a) ( ) 2x2xx2xxp 234 −++−= tem raízes i±± 1 ,1 . ( ) ( ) ( )1x2x2xxp 22 −⋅+−= Ex.: b) ( ) 10167 23 −+−= xxxxp tem raízes i±3 ,1 . ( ) ( ) ( )1x10x6xxp 2 −⋅+−= 5)Se ( )xp é de grau ímpar, então ( )xp possui ao menos uma raiz real. 6) Uma raiz 0x de ( )xp tem multiplicidade m se ( ) ( ) ( ) ( ) 0010"0'0 ===== − xpxpxpxp mK e ( ) 00 ≠xpm Ex.:1) 20 =x é raiz de multiplicidade 2 ( ) 8x6x2xp 23 +−= 2) 20 =x é raiz de multiplicidade 3 de ( ) 8465 234 −++−= xxxxxp 10 2 - Resolução de Equações Algébricas e Transcendentes 7) Valor numérico de um polinômio: para calcular, de forma usual, ( )ixp são necessárias n adições e ( ) 2 1+nn multiplicações. O Método de Horner faz esse cálculo com n adições e n multiplicações: ( ) 0121 parênteses 1 )))(((( axaxaxaxaxp nn n +++++= − − K321 KK Ex.: ( ) 4x)3x)1x)2x)4x3((((4x3xx2x4x3xp 2345 +−−+=−+−−+= 8) Se ( )xp é de grau n , então existe único polinômio de grau ( )xqn ,1− , tal que ( ) ( ) ( ) ( )αα pxqxxp +⋅−= . Se α é raiz de ( )xp então ( ) 0=αp . É o algoritmo de Briot-Ruffini utilizado para Deflacionar Raízes. Ex. ( ) ( )( ) ( ) ( )( ) ( ) ( )( ) ( ) 1483pe14846x10x3x3 22pe26x5x2x2 01pe10x6x1x1 10x16x7xxp 2 3 2 2 2 1 23 −=−−+−+⇒−= =++−−⇒= =+−−⇒= −+−= α α α 2.1.1-Enumeração das Raízes Enumerar as raízes de ( )xp é dizer quantas são as raízes e se positivas, negativas ou complexas. Regra de Descartes ou Regras de Sinais “O número de raízes reais positivas de ( )xp é igual ao número de variações de sinais na seqüência dos coeficientes ou menor do que este por um número inteiro par, sendo uma raiz de multiplicidade “m” contada como “m” raízes e não sendo considerados os coeficientes nulos”. O número de raízes reais negativas é obtido aplicando a regra de Descartes a ( )xp − Regra de Huat Se ( ) 00p ≠ e para algum k, 112 +− ×≤ kkk aaa então ( )xp possui raízes complexas. Regra da Lacuna Se ( ) 00 ≠p e para algum K, 0 e 0 11 >×= +− kkk aaa , então ( )xp tem raízes complexas. 11 E.B.Hauser – Cálculo Numérico Ex.: ( ) 153 345 −−−+= xxxxxp • ( )xp pode ter: 1 raiz real positiva, 2 raízes negativas e 2 complexas; ou • 1 raiz real positiva, nenhuma negativa, e 4 complexas. 2.1.2-Localização das raízes de p(x) Localizar as raízes de ( )xp é determinar a região do plano que contém todas as raízes. Cota de Cauchy: Toda raiz α (real ou complexa) de ( )xp satisfaz βα ≤ . Onde 0,lim 0 == ∞→ xxkkβ e n k n n k n nn k n n k a ax a ax a ax a ax 0122111 ++++= −−−−+ K Ex.: ( ) 13136,068,137,33 234 −+−+−= xxxxxp 4 1 23 1 )3136,068,137,33( +++=+ kkkk xxxx 00=x e 4 9575796715,340048339019,3 9552764745,336165556530,2 59519059204,361069395791,2 09469665820,364735839615,1 69397307712,3557483314773,0 195 184 173 162 151 ≤∴ == == == == == α MM xx xx xx xx xx 12 2 - Resolução de Equações Algébricas e Transcendentes 2.2- Separação de Raízes Reais de f(x)=0 a) Métodos Gráficos: Utiliza-se um dos seguintes processos: i) esboçar gráfico da função ( )xf e localizar as abcissas dos pontos onde a curva intercepta o eixo dos x. ii) de ( ) 0=xf obter uma equação equivalente ( ) ( )xfxf 21 = . Localizar no mesmo eixo cartesiano os pontos r onde as duas curvas se interceptem: ( ) ( ) ( ) ( ) ( ) 002121 =⇒=−⇒= rfrfrfrfrf b) Método Analítico: Seja ( )xf continua no intervalo [ ]ba, . Se ( ) ( ) 0<⋅ bfaf , então existe pelo menos uma raiz de f em ( )b,a . (Se o sinal de 'f é constante em ( )b,a a raiz é única nesse intervalo). Ex.: ( ) 393 +−= xxxp a) Análise gráfica: Logo, existem três raízes reais: ( ) ( ) ( )3,23r,1,02r3,41r ∈∈−−∈ b) analiticamente: Obs: Devemos dar uma atenção especial para os casos de: ¾ Raízes muito próximas. ¾ para raízes de multiplicidade par não ocorre troca de sinal. ( ) 373923,753113923,13325 323101334 −−−− −−−− xp x 13 E.B.Hauser – Cálculo Numérico Ex1: ( ) 3136,068,137,33 234 +−+−= xxxxxp 8,0 e 7,0 21 ≅= rr são raízes de multiplicidade 2. 2.3- Métodos para Resolução de equações algébricas e transcendentes Qualquer método deve observar um critério de parada, ao qual está associado um estimador de exatidão. Por exemplo, para onde LC ,,, 21 εε são dados: • C)x,x(DIGSE 1kk ≥+ • ( ) 1kxf ε≤ • 2 1 1 ε≤− + + k kk x xx • Lk > (número máximo de iterações) 2.3.1-Método da Bisseção ou Dicotomia (Algoritmo de quebra) Seja [ ] ℜ→baf ,: continua e tal que ( ) ( ) 0<⋅ bfaf . 1) Calcula-se o ponto médio 2 baxm += , dividindo-se [ ]ba, em dois novos intervalos : [ ]mxa, , [ ]bxm , 2) Se ( ) 0≠mxf e: i) ( ) ( ) 0<⋅ mxfaf então a raiz está em ( )mx,a . Volta-se para (1) ii) ( ) ( ) 0<⋅ bfxf m então a raiz está em ( )b,mx . Volta-se para (1) 3) Repete-se o processo até obter uma aproximação “razoável” da raiz, isto é, até que um critério de parada seja satisfeito. Características: É simples a convergência lenta mas garantida. A velocidade de convergência é DIGSE⋅3,0 /passo, isto é, a cada 3 ou 4 passos ganha-se um DIGSE . Ex: ( ) 11205,72 234 −−−+= xxxxxp a) Enumeração das raízes de ( )xp 4211 4031 total⊂ℜℜ −+ Regras de Huat e Lacuna não aplicam 14 2 - Resolução de Equações Algébricas e Transcendentes 1r 2r 3r 4r b) Cota de Cauchy: 0,11205,72 04 23 1 =+++=+ xxxx kkkk 46483007996,4 46481382609,481469532045,4 96478405782,467425618494,3 16472953921,460308011161,3 48211602868,1 17 164 153 142 1 = == == == = x xx xx xx x M M c) Separação de raízes ( ) 5,57617335495,3511515,8775,276 54321012345 −−−−− −−−−− xp x d) Cálculo da raiz ( )4,34r ∈ utilizando o método da bissecção. ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) M002874148,00352783203,313 20354003906,3;03515625,320152911151,020354003906,312 50356445312,3;03515625,30401045242,050356445312,311 0361328125,3;03515625,30897548754,00361328125,310 037109375,3;03515625,31891500395,0037109375,39 0390625,3;03515625,30095143543,003515625,38 039625,3;03125,33883184832,00390625,37 046875,3;03125,31900454168,1046875,36 0625,3;03125,3405335188,003125,35 0625,3;3871886353,20625,34 125,3;3660400391,9125,33 25,3;300390625,2525,32 35;39375,625,31 kIkxfkxk − − Obtemos 0354039062,312x4r =≈ com ( ) 096,413x,12xDIGSE ≅ 2.3.2-Método da Iteração Linear Consiste em escrever a equação ( ) 0=xf na forma ( )xGx = . Os pontos x* que satisfazem a condição ( )** xGx = são ditos pontos fixos de ( )xG e geometricamente representam os pontos de intersecção da reta xy = com a curva ( )xGy = . ( )xG é dita função de iteração do método. Inicia-se a iteração com um valor 0x próximo da raiz, e as outras aproximações são dadas por: ( ) K,2,1,0,1 ==+ ixGx ii 15 E.B.Hauser – Cálculo Numérico A seqüência de aproximação ix , converge para a solução x* da equação ( ) 0=xf sob certas condições . A construção de G não é única. A escolha de uma G apropriada é dita “problema do ponto fixo. Ex. 062 =−+ xx . ( ) 16, 1 6,6),6)6) 5432 2 1 −=+=−−=−=−= xGxGxGcxGbxxGa Embora não seja preciso usar métodos numéricos para encontrar as duas raízes reais 2,3 21 =−= αα e da presente equação, nota-se que: i) Tomando 1G e 5,10 =x , a seqüência }{ ix não converge para 2. ( )ii xGx 11 =+ ( ) ( ) ( ) ( ) ( ) M 8342,12078822 46095276,3475 00390625,59)0625,8(6 0625,8)75,3(6 75,3)5,1(6 415 314 2 213 2 112 2 011 −== −== −=−−== −=−== =−== xGx xGx xGx xGx xGx ii) Tomando 2G e 5,10 =x , a seqüência }{ ix converge para 2. ( ) ( ) ( ) ( ) ( ) ( ) ( ) M 10000298018,2x 69998807918,1x 50004768183,2x 39980924992,1x 40076263645,2x 9694363804,1x 61213203435,25,16x 627 526 425 324 223 122 021 == == == == == == =−== xG xG xG xG xG xG xG Teorema da Convergência: Seja α uma raiz isolada de f em [ ]ba, . Se i) G e G’são contínuas em [ ]ba, ; ii) ( )b,ax,1)x('G ∈∀〈 ; iii) ( ) ,...2,1,0k,b,a)kx(G1kxe0x =∈=+∈Ι , então a seqüência { kx }, gerada por )(1 kxGkx =+ ,converge para α. 16 2 - Resolução de Equações Algébricas e Transcendentes Ex: Utilizando o método da iteração linear calcule a raiz de ( ) 3xexf x += , com ( ) 5, 1 ≥+kk xxDIGSE ( ) 33 33 00 x x xx eex exxexf −=−=⇒ −=⇒=+⇒= Seja ( ) ( ) ( ) ( ) 33,00' 24,01'; 3 1' 3 3 ≅ ≅−−= −= G GexG exG x x * G e G’ são continuas em [-1,0] e ( ) ]0,1[1' −∈∀< xxG . Logo , a seqüência gerada por 31 ix i ex −=+ converge para α ]0,1[−∈∀x . Seja 5,00 −=x M 772882595,0773204044,0 772884374,0771636903,0 772877469,0777723518,0 772904269,0754152577,0 772800243,0846481725,0 105 94 83 72 61 −=−= −=−= −=−= −=−= −=−= xx xx xx xx xx ( ) 34,5),( 000003188,0 109 9 ≅ = xxDIGSE exf *G não tem Maximo nem mínimo local em [0,1], testa-se então só os extremos. Características do Método da Iteração Linear: ¾ Não garante a convergência para toda função continua. ¾ Necessita do calculo de G’(x). ¾ Pode ocorrer dificuldade para encontrar G(x). ¾ A convergência é linear para raízes simples (a cada passo do método o erro é reduzido por um fato constante). ¾ A velocidade de convergência depende de ( )xG' , quanto menor este valor, maior será a convergência. 17 E.B.Hauser – Cálculo Numérico Obs.: Dada a equação ( ) 0=xf , existem infinitas funções G(x) que são funções de iteração. A forma geral destas funções é: ( ) )()( xfxAxxG ⋅+= , com a condição que em x*, ponto fixo de G(x), se tenha 0)( * ≠xA . Temos que: *** )(0)( xxGxf =⇔= Com efeito: (⇒ ) Seja *x tal que 0)( * =xf . ******* 0)()()()( xxAxxfxAxxG =⋅+=⋅+= 0)(0)(0)()( )()()()( **** ****** ≠=⇒=⋅⇒ =⋅+⇒=⇐ xApoisxfxfxA xxfxAxxxGSe 2.3.3-Métodode Newton-Raphson Procura garantir e acelerar a convergência do método da Iteração Linear, escolhendo para a função de iteração a função G(x) tal que G’(x)=0 Dada a função 0)( =xf , tomamos a forma geral para G(x): ),(')(1)(')(')()()('1)(' )(')()()('1)(')()()( ******** xfxAxGxfxAxfxAxG xfxAxfxAxGxfxAxxG ⋅+=⇒⋅+⋅+= ⇒⋅+⋅+=⇒⋅+= pois x * é ponto fixo de G )0)()(( *** =⇒= xfxGx . Agora , )(' 1)(0)(')(10)( * **** xf xAxfxAxG −=⇒=⋅+⇒= . Tomemos )(' )()( )(' 1)( xf xfxxGe xf xA −=−= . Então dada ( )xf , ( )xG é tal que ( ) 0' * =xG , pois ( ) 22 2 )]('[ )(")(' )]('[ )(")()]('[1' xf xfxf xf xfxfxfxG ⋅=⋅−−= e .0)('0)('0)( *** ≠=⇒= xfsexGxf Portanto, iniciando-se a iteração com um valor 0x escolhido, a seqüência )( kx é determinada por: 0)(',...2,1,0, )(' )( 1 ≠=−=+ kxfk kxf kxfkxkx Geometricamente , conforme podemos observar na próxima figura, dado kx , 1+kx é abscissa do ponto de intersecção da reta tangente à curva ( )xf em ))(,( kk xfx e o eixo dos “x”. Assim, )(' )( )(' )( 1 1 k k kkk kk k xf xf xxxf xx xf tg −=⇒=−= ++ θ 18 2 - Resolução de Equações Algébricas e Transcendentes Convergência: (é trabalhoso mostrar que 1)x('G < ). O método de Newton-Raphson converge se: 2 2 ))('()(")(1))('( )(")()(' xfxfxf xf xfxfxG <⇒<= . Para raízes simples a convergência é quadrática e para raízes duplas ou triplas é linear. Escolha do ponto inicial: Seja )b,a(∈α raiz de f . Se bxbfbf axafaf =⇒>⋅ =⇒>⋅ 0 0 0)(")( 0)(")( Caso contrário, pode-se considerar 2 )( 0 bax += . Ex.: 1) Estimar o valor da única raiz real de xlnx2)x(f += , utilizando Newton-Raphson. 1kx + θ kx 1kx + )x(f 1k+ )x(f k xlnx2)x(f += x 12 xlnx2x,x − +− 19 E.B.Hauser – Cálculo Numérico k kk k1k x 12 xlnx2xx + +−=+ 34 3 2 1 0 xx 426302751,0x 42699599,0x 42,0x 5,0x = = = = = Logo a aproximação para a raiz é 426302751,0=α com 9 dígitos significativos corretos. 2) Calcular a raiz ]4,3[4 ∈r do polinômio dado anteriormente: 11205,72)( 234 +−+= xxxxxp . 40)4(")4(0)3(")3( 0 =⇒>⋅<⋅ xppepp 201564 11205,72 23 234 1 −−+ −−−+−=+ kkk kkkk kk xxx xxxxxx 56 5 4 3 2 1 0 03524990,3 03525211,3 03709653,3 08982331,3 36397059,3 4 xx x x x x x x = = = = = = = Obs: O Método de Newton pode divergir devido ao uso da tangente, oscilando indefinidamente. Uma aproximação para a raiz é 03524990,34 =r com 9 dígitos significativos e 5 iterações. 20 2 - Resolução de Equações Algébricas e Transcendentes Isto acontece quando: i) Não há raiz real ii) Ocorre simetria de ( )xf em torno de α iii) O valor inicial 0x está longe da raiz exata, fazendo que outra parte da função prenda a iteração. 2.3.4-Método da Secante Uma desvantagem do Método de Newton-Raphson é o calculo do valor numérico de ( )xf ' a cada iteração. O método da secante contorna este problema, substituindo a derivada pelo quociente das diferenças: 1 1 )()()(' − − − −≅ kk kk k xx xfxf xf onde kx e 1−kx são duas aproximações para a raiz de ( )xf . A formula iterativa é: )()( )()( 1 1 1 − − + − ⋅−= kk kkk k xfxf xfxx x Geometricamente, a partir das aproximações para a raiz de kx e 1+kx , o ponto 1+kx é dado pela abscissa do ponto de intersecção do eixo Ox e da reta secante que passa por ))(,())(,( 11 kkkk xfxexfx −− . Características do método da secante: 9 Garante a convergência para toda função continua 9 Pode divergir se )()( 1−≅ kk xfxf 9 Convergência mais lenta que o Método de Newton e mais rápida que Bissecção e Iteração Linear. 9 São necessárias duas aproximações da raiz a cada iteração. Ex: )0,1(e21x172x53x)x(p −∈++−= α )21175()21175( )21175()( 1 2 1 3 1 23 23 1 1 ++−⋅++− ++−⋅−−= −−− − + kkkkkk kkkkk kk xxxxxx xxxxxxx 67 6 5 4 3 2 1 10 629321148566,0 89321148568,0 069321123947,0 869317876515,0 599601423487,0 898888888888,0 1,1 xx x x x x x x xx = −= −= −= −= −= −= =−= 21 E.B.Hauser – Cálculo Numérico Exercícios 1) Uma partícula é arremessada verticalmente, a partir do solo, com uma velocidade inicial ov .Desprezando a resistência do ar, supomos que a posição p partícula é dada por: 2 2 tgtov)t(d −= , onde g é aceleração da gravidade. Determinar a altura máxima atingida pela partícula e o instante em que ocorreu. 2) Uma corrente oscilante num circuito elétrico é descrita por I = )t2sen(te9 π− , t em segundos. Determinar todos os valores positivos de t para os quais I = 3.5. 3) A concentração de uma bactéria poluente num lago é descrita por C = t075.0e5.2t5.1e70 −+− Encontrar o tempo para que a concentração seja reduzida para nove. 4) O deslocamento de uma estrutura está definido pela seguinte equação D = wtcoskte8 − onde k = 0.5 e w = 3. a) Determinar graficamente uma estimativa inicial do tempo necessário para o deslocamento decrescer para 4. b) Usar o método de Newton-raphson para determinar essa raiz 5) Enumerar, localizar, separar e calcular (via Newton-Raphson e/ou Bissecção ), se possível, todas as raízes dos polinomios tendo como critério de parada DIGSE (xk , xk+1) ≥ 5. No caso de raízes múltiplas, determinar a multiplicidade da raiz e calcular as demais utilizando deflação. a) p(x) = x x x3 22 3 5− + − b) p(x) = x5 - 15,5x4 +77,5x3 - 155x2 +124x -31 c) p(x) = x x x x4 3 2121 2247 15043 34300− + − + d) p(x) = x x x x4 3 2115 1575 7625 12500− + − + e) p(x) = x x x x4 3 23 337 168 0 3136− + − +. . . f) p(x) = x4-11,101x3+11,1111x2-1,011x+0,001 g) p(x) = x9- 4x8 + 3,9x7 +3,02x6 - 5,5295x5 - 0,84732x4 + 2,83536x3 + 0,37152x2 -0,59616x - 0,15552 h) p(x) = x3 - 2081,93x2 +1424,64x- 3,19728 i) p(x) = x4 + 1,98x3 +1,1424x2 +0,162602x - 0,00174225 6) Localizar graficamente e calcular ( via Newton-Raphson e/ou Método da Iteração Linear) todas as raízes, com DIGSE(xk , xk+1) ≥ 5: 22 2 - Resolução de Equações Algébricas e Transcendentes a) f(x) = x2 + ln x b) f(x) = x + 2 cos x c) f(x) = 2x + ln x d) f(x) = cos x + ln x + x e) f(x)= x + e Bx− 2 para B = 1,5,10,25,50 7) Responder resumidamente: a) Em que consiste o princípio da bissecção para determinar a primeira aproximação de uma raiz de uma equação f(x)=0? b) Explicar o método da iteração linear para calcular uma raiz de uma equação f(x)=0, partindo de uma primeira aproximação x0. c) Quando não converge a iteração linear? d) Quando não converge o Método de Newton Raphson? e) Interpretar geometricamente o Método de Newton-Raphson. f) Qual a vantagem de se utilizar o Algoritmo de Horner para se avaliar o valor do polinomio num ponto? Exemplificar. Respostas: 1) O deslocamento máximo é vo2/2g e ocorreu em t = vo/g. 2) t = 0.06835432097 e t = 0.4013436927 3) t = 1.556787935 23 E.B.Hauser – Cálculo Numérico 4) t = 0 .3151660803 5-a) r+ r- ¢ T 3 0 0 3 1 0 2 3 p tem raízes complexas. Existe uma raiz real em (1,2) Raízes:x=1,84373427779 x= 0.07813286110 ± 1.644926378i . b) Raízes: .4555300547, 1.092601944, 1.940878206, 4.011783883, 7.999205912 c) Raízes: R1=100 e R2= 7(multiplicidade 3), não tem raízes complexas. d) Raízes: R1=10 e R2=5(multiplicidade 3), não tem raízes complexas. e) Raízes: R1=0.7(multiplicidade 2), R2=0.8(multiplicidade 2) f) Raízes: R1=0,001 R2=0,1 R3=1 R4=10 g) Raízes: R1=-0,5(multiplicidade 4), R2=1,2(multiplicidade 5) h) Raízes: R1= 0.002251681490, R2 = 0.6822607762 , R3= 2081.245488 i) Raízes: R1=-1,01 R2=-0,75 R3=-0,23 R4=0,01 6-a) R≅ .6529186400 b) R≅ -1.029866529 c) R≅ 0 .42630275100 d) R≅ .2875182754 e) Existe única raiz de f em (-1,0) 24 3. Sistemas de Equações Lineares O sistema com n equações lineares e n variáveis bxaxaxaxa bxaxaxaxa bxaxaxaxa nnnn33n22n11n 2nn2323222221 1nn1313212111 =++++ =++++ =++++ L MMMMM L L pode ser representado matricialmente por BAX = , onde ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎣ ⎡ = aaa aaa aaa A nn2n1n n22221 n11211 K MMM K K , ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎣ ⎡ = x x x X n 2 1 M , ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎣ ⎡ = b b b B n 2 1 M e A é a matriz dos coeficientes, X é o vetor das incógnitas e B o vetor dos termos independentes. 3.1- Introdução à problemática de sistemas Um SEL pode ser mal condicionado, bem condicionado ou sem solução. Um sistema é “mal condicionado” se uma pequena alteração nos dados pode provocar uma grande alteração na solução. Por exemplo: (a) ⎩⎨ ⎧ =+ =+ 5yx 95,4y98,0x tem solução exata ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 5,2 5,2 x (b) ⎩⎨ ⎧ =+ =+ 5yx 95,4y99,0x tem solução exata ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 0,5 0,0 x Uma alteração de 1% (0,01 no coeficiente 0,98) ocasionou uma alteração de 100% na solução. No caso de um sistema linear de ordem 2, cada equação representa uma reta. Resolver o sistema significa determinar a intersecção das duas retas. Três casos são possíveis: Bem condicionado Não tem solução. Mal condicionado 0det ≠ 0det = 0det ≅ retas concorrentes retas paralelas (perturbação em 2ℜ ) R1 R2 R2 R1 R2 R1 25 E.B.Hauser – Cálculo Numérico Exemplo2: (a) ⎩⎨ ⎧ =+ =+ 503,25y501,7x5,1 17y5x tem solução exata: ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 3 2 x (b) ⎩⎨ ⎧ =+ =+ 501,25y501,7x5,1 17y5x tem solução: ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 1 12 x 3.2-Medidas de Condicionamento O determinante normalizado da matriz dos coeficientes A é dado por 2 kn 2 2k 2 1kk n21 aaaondeAdetANORM LL ++== αααα , para k = 1, 2, ..., n ( )1,1ANorm −∈ e quanto mais afastado de 1± (isto é, quanto mais próximo de zero) mais mal condicionado é a matriz A. Retomando o Ex2: 2) ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 501,75,1 51 A ( ) ( ) 586495509853,7501,75,1 90990195135,5251 2 1 22 2 2 1 1 =+= =+= α α 38740000256377,0 0050000128,39 001,0001,0Anorm 001,05,15501,7Adet 21 ==⋅= =⋅−= αα Agora, como pode ser medido o condicionamento de um sistema linear? Dado um SEL BAX = , o seu número de condicionamento é dado por: 1AA)A(Cond −= . Quanto maior for )A(Cond , mais sensível será o sistema linear. Utilizamos ∑== =≤≤∞ n 1i ij ni1 amaxAA , a norma do máximo das linhas. Ex: ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛= 501,75,1 51 A , ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − −=− 10001500 50007501 A 1 51 1 101501,112521AA)A(Cond 12501A,001,9A ⋅≅== == − − 26 Sistemas de Equaçõe Lineares 3.3-Método de Resolução de Sistemas Métodos Diretos: A solução exata é obtida realizando-se um número finito de operações aritméticas em ℜ (com precisão infinita): Eliminação de Gauss e Gauss Jordan. Métodos Iterativos: A solução x é obtida como limite de uma seqüência de aproximações sucessivas x1, x2, ... . Método de Eliminação de Gauss Algoritmo básico de Gauss: A solução de BAX = é calculada em duas etapas: 1º- Triangularização : Mediante operações elementares nas linhas, a matriz A é transformada numa matriz triangular superior. Algoritmo: para 1n,,1k −= K (indica a linha do pivô) para n,,1ki K+= (indica a linha a transformar de A) kk ik a am −= 0aik = para nkj ,,1K+= (indica a coluna a transformar da linha i) amaa kjijij ⋅+= bmbb kii ⋅−= Obs.: Se 0aii = é necessário trocar de linha, se possível. 2º- Retro-substituição: A triangularização transforma o sistema BAX = , no sistema equivalente: dxc dxcxc dxcxcxc dxcxcxcxc nnnn 3nn3333 2nn2323322 1nn1313212111 = =++ =+++ =++++ LLLLLLL L L L cuja solução é dada por: ( ) ( ) 11 nn13132121 1 1n,1n nn,1n1n 1n nn n n a xcxcxcdx , a xcd x, c dx +++−= −== −− −−− L Teorema: O método de Gauss produz sempre a solução exata do sistema BAX = (utilizando precisão infinita) se 0det ≠A e as linhas de A forem permutadas quando 0aii = . 27 E.B.Hauser – Cálculo Numérico Quando é utilizada aritmética do ponto flutuante, erros de arredondamento podem comprometer a solução obtida. Ex.: Em ( )98,98,3,10F − , com arredondamento para número mais próximo de máquina “ox”, a elminação de Gauss aplicada ao sistema ⎩⎨ ⎧ =+ =+ 6y2x2 5y2x0002,0 produz 0x = e 5,2y = o que não satisfaz a segunda equação do sistema. (Obs.: multiplicador 000.10m −= por ( )mLLL 122 += ) Gauss com Pivotamento Parcial O método consiste em trocar linhas (ou colunas) de maneira a minimizar a propagação de erros nas operações. Escolhas dos pivôs: 1º pivô é o elemento de maior valor absoluto da coluna 1. 2º pivô é o elemento de maior valor absoluto da coluna 2 e da diagonal para baixo. Procede-se da mesma forma para os demais pivôs: Aplicando a técnica ao último exemplo ⎩⎨ ⎧ =+ =+ 5y2x0002,0 6y2x2 , em ( )98,98,3,10F − , com arredondamento para número mais próximo de máquina “ox”, obtemo 5,0x = e 5,2y = . (Obs.: Neste caso o multiplicador é menor m= - 0,0001) Método de Gauss-Jordan (Matriz Inversa) A solução do SEL BAX = é calculado utilizando BAX 1−= se 0Adet ≠ . Obs.: Ver exercício 9. Métodos Iterativos Os sistemas lineares de grande parte são, em geral, esparsos (muito elementos 0aij = ). Os métodos diretos não preservam a esparsidade, enquanto que os métodos iterativos sim, além de apresentarem relativa insensibilidade ao crescimento dos erros de arredondamento. No sistema original BAX = , supondo n,,1i,0aii K=≠ , o vetor X é isolado mediante a separação de diagonal principal: 1º pivô 2º pivô 3º pivô 28 Sistemas de Equaçõe Lineares ( ) ( ) ( )1n1n,n22n11nn nn n nn23231212 22 2 nn13132121 11 1 xaxaxab a 1x xaxaxab a 1x xaxaxab a 1x −−−−−−= −−−−= −−−−= K MM K K Metodo de Gauss-Jacobi Dada a aproximação inicial X0, o processo iterativo produz aproximações sucessivas KK ,,,, 21 kXXX , obtidas de: ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ( ) ⎟⎠⎞⎜⎝⎛ −−−−= −−−−= ⎟⎠ ⎞⎜⎝ ⎛ −−−−= −−+ + + k 1n1n,n k 22n k 11nn nn 1k n k nn2 k 323 k 1212 22 1k 2 k nn1 k 313 k 2121 11 1k 1 xaxaxab a 1x xaxaxab a 1x xaxaxab a 1x K MM K K Métodode Gauss-Seidel Para X0 dado: ( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ( )( )1k 1n1n,n1k22n1k11nn nn 1k n k nn3 )k( 434 1k 232 1k 1313 33 1k 3 k nn2 k 323 1k 1212 22 1k 2 k nn1 k 313 k 2121 11 1k 1 xaxaxab a 1x xaxaxaxab a 1x xaxaxab a 1x xaxaxab a 1x + −− +++ +++ ++ + −−−−= −−−−−= −−−−= −−−−= K MM K K K Critério de Convergência O teorema abaixo estabelece uma condição suficiente para a convergência dos métodos de Gauss-Jacobi e de Gauss-Seidel. Teorema Dado o sistema linear YAX = , se a matriz A é Diagonalmente Dominante, isto é, se iaa ij ijii ∀∑> ≠ , então tanto o método de Jacobi como o de Gauss-Seidel gera uma seqüência ( )( )kX convergente para a solução do sistema dado, independente da escolha da aproximação inicial ( )0X . 29 E.B.Hauser – Cálculo Numérico Exemplo: Resolver o sistema abaixo por Gauss- Jacobi e Gauss-Seidel. Critério de Parada: erro absoluto da solução calculada for menor que 10-3. ⎪⎪⎩ ⎪⎪⎨ ⎧ − − =+ =+ =+− =−+ 5,1 1 3 2 xx5,0 xx2 x2x xx4x 31 41 43 421 Reordenamos a fim de satisfazer ao critério de convergência. ⎪⎪⎩ ⎪⎪⎨ ⎧ − − =+− =+ =−+ =+ 3 5,1 2 1 x2x xx5,0 xx4x xx2 43 31 421 41 12 5,01 114 12 −> > −+> > Como a matriz dos coeficientes , após a reordenação, é diagonalmente dominante, então a aplicação dos métodos de Gauss-Jacobi e Gauss-Seidel produzirá uma seqüência ( )( )kX convergente para a solução exata. Gauss-Jacobi Fórmulas iterativas: ( ) ( ) ( ) ( ) ( ) ( ) ( )( ) ( ) ( )k 3 )k( 31k 4 k 1 1k 3 k 4 k 1 )k( 4 )k( 11k 2 k 4 )k( 41k 1 x5,05,1 2 )x3( x x5,05,1x x25,0x25,05,0 4 )xx2( x x5,05,0 2 x1 x +−=+−= −= +−−=+−−= −= ⎟⎠ ⎞⎜⎝ ⎛ − = + + + + Aproximação inicial: ( ) 0X 0 = . 30 Sistemas de Equaçõe Lineares k x1 x2 x3 x4 0 0 0 0 0 1 0,5 0,5 1,5 -1,5 2 1,25 -1,0 1,25 -0,75 3 0,875 -1,0 0,875 -0,875 4 0,9375 -0,9375 1,0625 -10625 5 1,03125 -1,0000 1,03125 -0,96875 6 0,984375 -1,0000 0,984375 -0,984375 7 0,9921875 -0,9921875 1,0078125 -1,0078125 8 1,00390625 -1,0000 1,00390625 0,9960375 9 0,998046875 -1,0000 0,998046875 -0,998046875 10 0,9990234375 -0,999023437 1,000976563 -1,000976563 11 1.000488282 -1.0000000 1.000488281 -0,999511718 : : : : : 40 1 -1 1 -1 na 12a. iteração consegue-se ( ) 51112 10xx −<− Gauss-Seidel Fórmulas iterativas: ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )1k 3 1k 4 1k 1 1k 3 k 4 1k 1 1k 2 k 4 1k 1 x5,05,1x x5,05,1x x25,0x25,05,0x x5,05,0x ++ ++ ++ + +−= −= +−−= −= k x1 x2 x3 x4 0 0 0 0 0 1 0,5 -0,625 1,25 -0,875 2 0,9375 0,953125 1,03125 -0,984375 3 0,9921875 -0,994140625 1,00390625 0,99806875 4 0,9990234375 -0,9992675781 1,000488281 -0,999755895 5 0,999877929 -0,999908447 1,000061035 -0,999969482 : : : : : 12 1 -1 1 -1 ( ) ( ) 545 10xx −<− 31 E.B.Hauser – Cálculo Numérico EXERCÍCIOS 1) Uma consideração importante no estudo da transferência de calor é a de se determinar a distribuição de temperatura numa placa, quando a temperatura nas bordas é conhecida. Supomos que a placa da figura represente um seção transversal de uma barra de metal, com fluxo de calor desprezível na direção perpendicular à placa. Sejam T1, T2 , ..., T6 as temperaturas nos seis vértices interiores do reticulado da figura. A temperatura num vértice é aproximadamente igual à média dos quatro vértices vizinhos mais próximos (à esquerda, acima, à direita e abaixo). Por exemplo: T1 = ( 10 + 20 + T2 + T4 )/4 ou 4T1 = 10 + 20 + T2 + T4 a) Escrever o sistema de seis equações, cuja solução fornece estimativas para as temperaturas T1, T2 , ..., T6. b) Resolver o sistema utilizando o sistema MapleV. 2) Num experimento num túnel de vento, a força sobre um projétil devido à resistência do ar foi medida para velocidades diferentes. 119 10 3.74 8 6.39 6 8.14 4 90.2 2 0 0 força velocidade Expressar a força como função da velocidade aproximando-a a um polinômio de quinto grau: 5 5 4 4 3 3 2 21o vavavavavaa)v(f +++++= Verificar a validade de )(vf encontrada e obter uma estimativa para a força sobre o projétil quando ele estiver se deslocando a uma velocidade de 1, 3, 5, 7 e 9 unidades de velocidade. 3) Considere o sistema (#) ⎪⎩ ⎪⎨ ⎧ =++ =++ =++ 78,0x2,0x25,0x33,0 08,1x25,0x33,0x5,0 83,1x33,0x5,0x 321 321 321 . a) (#) é bem ou mal condicionado? Porque? O que isso significa? b) Resolver (#) pelo método de Gauss sem pivotamento. 1 20o 20o 20o 20o 20o 20o 40o 40o 10o 10o 1 2 3 4 5 6 32 Sistemas de Equaçõe Lineares 4) Idem ao 3 para ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ =++++ =++++ =++++ =++++ =++++ 2520/1879x9/1x8/1x7/1x6/1x5/1 840/743x8/1x7/1x6/1x5/1x4/1 420/459x7/1x6/1x5/1x4/1x3/1 60/87x6/1x5/1x4/1x3/1x2/1 60/137x5/1x4/1x3/1x2/1x 54321 54321 54321 54321 54321 5) Resolver por Eliminação de Gauss com pivotamento parcial. a) ⎪⎪ ⎪⎪ ⎩ ⎪⎪ ⎪⎪ ⎨ ⎧ =+−++−+ −=−−−+−+ =+++−++ =−−++++− =++−++− =+−++++ =+−−++− 34x9xx7x2xx2x7 10xx4xx11x8x9x2 13xxx2x7xx3x4 8x2xx7x2x3xx 24x11x6xx5x7xx2 17x2x9xx15x9xx 47x8x3x2x4x6x5x3 7654321 7654321 7654321 7654321 7654321 7654321 7654321 b) ⎪⎩ ⎪⎨ ⎧ =−−+ =−− =++− 255x10x35x 4414x41x72x 8342x2xx 321 321 321 c) ⎪⎪⎩ ⎪⎪⎨ ⎧ =−++ −=++− =−++ =+−+ 9.9xx5.2x2.0x 9.3x2.5xxx3.0 9.21x5.8x4x5.0x4.0 7.2xx1.0xx2 4321 4321 4321 4321 6) Referente ao sistema linear AX=B, verificar se a afirmação é Verdadeira ou falsa . i) Se det A = 0 então o sistema não tem solução .Justificar verificando o que acontece em : a) ⎪⎩ ⎪⎨ ⎧ =+− =+ =+− 2xx2 4xx 5xxx 32 21 321 e b) ⎪⎩ ⎪⎨ ⎧ =+ =− =+ 0xx 0xx 0xx 32 31 21 ii) Se A não é uma matriz Diagonal Dominante então os métodos de Gauss-Jacobi e Gauss-Seidel não geram uma sequência convergente para a solução . Justificar verificando o que acontece com a) ⎩⎨ ⎧ −=− =+ 3y3x 3yx b) ⎩⎨ ⎧ =+ −=− 3yx 3y3x 7) Resolver pelo Método de Gauss-Seidel. Apresentar as fórmulas iterativas e uma garantia de convergência (se possível). a) ⎪⎪⎩ ⎪⎪⎨ ⎧ =+−+ =−++ =+++ =++− 82x7x12x66x17 34x18x11x22x56 26x10x17x16x11 18x20x47x5x3 4321 4321 4321 4321 b) ⎪⎩ ⎪⎨ ⎧ =++ =+− =++ 31x6x3x2 1x2x9x 30x3x2x8 321 321 321 33 E.B.Hauser – Cálculo Numérico 8) Resolver o sistema de equações algébricas não lineares: ⎪⎪⎩ ⎪⎪⎨ ⎧ =++ =++ =++ 15z4yx 18zy4x 11zyx4 2 2 2 9) Resolver ⎪⎩ ⎪⎨ ⎧ =++ =++ =++ 3321 2321 1321 bx4x3x5 bx2x3x bx7xx2 para: a) b1 =16 b2 = -5 b3=11 b) b1 =25 b2 = -11 b3 = -5 c) b1 =3 b2 = 5 b3 = -5 10) Utilizando Eliminação Gaussiana calcular detA. a) A = ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎣ ⎡ −−− −−− −− −−−− 1322 1413 2101 2132b) A = ⎥⎥ ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎣ ⎡ − 01301 25.09110 30350 00023 05.0807 11) Qual o Resíduo produzido pela solução aproximada X’= [ -3 4 0]T de ⎪⎩ ⎪⎨ ⎧ =++ =++ =++ 64.0x25.0x21.0x15.0 52.0x24.0x16.0x12.0 84.0x12.0x36.0x24.0 321 321 321 Respostas: 1) Solução obtida utilizando o MapleV: >solve({4*T1=10+20+T2+T4,4*T2=T1+20+T3+T5,4*T3=T2+20+40+T6, 4*T4=10+T1+T5+20,4*T5=T4+T2+T6+20, 4*T6=T5+T3+40+20}, {T1,T2,T3,T4,T5,T6}); {T6 = 190/7, T5 = 150/7, T4 = 120/7, T3 = 190/7, T1 = 120/7, T2 = 150/7} > evalf(%); {T6 = 27.14285714, T5 = 21.42857143, T4 = 17.14285714, T3 = 27.14285714, T1 = 17.14285714, T2 = 21.42857143} 2) Solução utilizando o sistema MapleV: >solve({a0=0,2.90=a0+a1*2+a2*4+a3*8+a4*16+a5*32, 14.8=a0+a1*4+a2*(4^2)+a3*(4^3)+a4*(4^4)+a5*(4^5), 39.6=a0+a1*6+a2*(6^2)+a3*(6^3)+a4*(6^4)+a5*(6^5), 74.3=a0+a1*8+a2*(8^2)+a3*(8^3)+a4*(8^4)+a5*(8^5), 119=a0+a1*10+a2*(10^2)+a3*(10^3)+a4*(10^4)+a5*(10^5)}, {a0,a1,a2,a3,a4,a5}); {a0 = 0, a2 = -1.194791667, a5 = .002604166667, a4 = -.07005208333, a3 = .6614583333, a1 = 1.712500000} >f:=x->1.712500000*x-1.194791667*x^2+.6614583333*x^3-.7005208333e- 1*x^4+.2604166667e-2*x^5; > validade:=[f(0),f(2),f(4),f(6),f(8),f(10)]; 34 Sistemas de Equaçõe Lineares validade := [0, 2.899999998, 14.80000000, 39.60000000, 74.29999994, 119.0000000] > estimativas:=[f(1),f(3),f(5),f(7),f(9)]; estimativas := [1.111718750, 7.202343750, 25.73046873, 55.89609367, 94.9992188] 3-a) Mal condicionado. NORM A ≅ 0,000181. Uma pequena perturbação nos dados de entrada pode causar uma grande alteração na solução. b) A solução exata é X=[1 1 1]T 4 - NORM A ≅ 0,00257, a solução exata é X=[1 1 1 1 1 1 1]T 5- a)X=[-21,86188 11,46568 2,376447 -8,514801 0,7475478 -15,50981 18,08498]T b) solução exata X=[1 0 2]T c) solução exata X=[1 2 3 -1]T 6-i) a) detA=0 e o sistema e incompatível b) detA = 0 e o sistema tem infinitas soluções dadas por x=z e y=-z ii) A não é matriz Diagonal Dominante e Gauss-Jacobi e Gauss-Seidel a) converge (oscilando) para a solução exata [1.5 1.5]T b) diverge 7 - a) X(10) = [-0.930569 1.901519 1.359500 -1.906078]T X(35) = [-1.076888 1.990028 1.474477 -1.906078]T b) solução exata X = [2 1 4]T 8 - A solução exata é x =1 y= 2 z = 3. 9 - soluções exatas: a) X = [3 -4 2]T b) X = [2 -7 4]T c) X = [-3 2 1]T 10 - a) detA = -55 b) det A = 706.875 11 - R = [0.12 0.24 0.25]T 35 4.Interpolação Polinomial 4.1 Objetivo: Dado um conjunto de n+1 pontos distintos ))x(f,x( ii , n,...,1,0i = , queremos determinar o polinômio p(x) talque ),x(f)x(p ii = e para os demais pontos do intervalo )x(f)x(p ≅ , [ ]no x;xx∈∀ . O polinômio p(x) é dito polinômio aproximante ou interpolador de f(x) no intervalo [ ]nxox , . 4.2 Aplicações • Obter uma expressão analítica aproximada de uma função que é conhecida em apenas um número finito de pontos; • Avaliar a função num ponto não tabelado [ ]nxoxx ;*∈ ; • Determinar aproximações para ∫ nxox dxxf )( , substituindo a f(x) pelo polinômio interpolador; • Calcular uma aproximação para )(' xf para [ ]nxoxx ;∈ , substituindo f(x) por p(x). 4.2 Existência e Unicidade da Solução Dados ℜ∈ix e ,)x(f i ℜ∈ n,...,1,0i = , procuramos nP)x(p ∈ tal que iii x),x(f)x(p ∀= . Seja ∑ = =++++= n k kxka nxnaxaxaoaxp 0 ...221)( . 36 E.B.Hauser – Cálculo Numérico Então de =)x(p i ∑ = = n 0k i k ik ),x(fxa obtemos: ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ =++++ =++++ =++++ n n nn 2 n2n10 1 n 1n 2 12110 0 n 0n 2 02010 fxa...xaxaa .............................................. fxa...xaxaa fxa...xaxaa , o qual representa um sistema linear de ordem n+1 onde as n+1 incógnitas são os ,ka n,0k = e a matriz dos coeficientes é dada por: A = ⎥⎥ ⎥⎥ ⎥⎥ ⎥⎥ ⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎢⎢ ⎢ ⎣ ⎡ n n 2 nn n 2 2 22 n 1 2 11 n 0 2 00 x..xx1 ...... ...... ...... x..xx1 x..xx1 x..xx1 De acordo com Vandermonde, det A = ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ −∏∏ += − = n 1ij ij 1n 0i )xx( . Como os pontos são distintos, a diferença ixjx − será sempre diferente de zero, e então detA ≠ 0. Portanto o polinômio interpolador existe e também é único. 4.3Polinômio Interpolador de Newton Para Diferenças Finitas Ascendentes Dados ( )ii y,x , )x(fy ii = i=0,1,2,...,n satisfazendo hxx....xxxx 1nn1201 =−==−=− − , isto é , hxx i1i =−+ . Para k= 1, 2, ...n, e i= 0, 1, 2, ..., n-k , a diferença finita de ordem k é dada por i 1k 1i 1k i k yyy −+− −= ∆∆∆ . E, o polinômio interpolador de Newton para diferenças finitas ascendentes é dado por : o n n 1n1o o 2 2 1o 0 o o y h!n )xx()xx)(xx( y h!2 )xx)(xx( y h )xx( y)x(p ∆∆∆ −−−−++−−+−+= LL 37 4-Interpolação Polinomial Exemplo1: O alongamento de uma mola foi medido em função da carga aplicada. Obteve-se: 3,6 8 0,5 6 5,2 4 0,1 2 )cm(oalongament )kg(aargc Estimar o alongamento para o caso de ser aplicada uma carga de 7kg , utilizando o polinômio interpolador de Newton para diferenças finitas. Solução: 1) Tabela das diferenças finitas: i ix iy iy∆ i2 y∆ i3 y∆ 0 2 1 1.5 1 -2.2 1 4 2.5 2.5 -1.2 -------------------- 2 6 5 1.3 -------------------- -------------------- 3 8 6.3 -------------------- -------------------- -------------------- 2) Polinômio interpolador: ( ) ( )( ) ( )( )( ) )2,2.( )2(!3 6x4x2x1. )2(!2 4x2x5,1. 2 2x1)x(p 32 − −−−+−−+−+= 7,2x016666667,2x675,0x30458333333,0)x(p 23 +−+−= 3) Verificação de validade de p(x) : 9999999994,0)2(p = ( 1≅ ) 4999999999,2)4(p = ( 5,2≅ ) 5)6(p = e p(8)=6.3 4) Estimativa do alongamento ao se aplicar uma carga de 7kg: O alongamento da mola neste caso é aproximadamente 9375,5)7(p = cm. 5) Análise gráfica do problema: 38 E.B.Hauser – Cálculo Numérico 4.4 Polinômio Interpolador de Newton Para Diferenças Finitas Divididas Dados ( )iyix , , ),( ixfiy = i= 0, 1, 2, ..., n os pontos ix podem ter um espaçamento qualquer, não necessariamente eqüidistantes. O polinômio de Newton para diferenças divididas é dado por: oy n nxxxxoxxoyxxoxxyoxxoyxp ∆−−−−++∆−−+∆−+= )1)...(1)((... 2)1)((0)()( , onde, para k = 1, 2, ..., n, e i= 0, 1, 2, ..., n-k a diferença dividida de ordem k é dada por iki i 1k 1i 1k i k xx yy y − − = + −+− ∆∆∆ Por exemplo, para o caso de n=4: i ix iy iy∆ i2 y∆ i3 y∆ i4 y∆ 0 0x 0y 01 01 xx yy − − 02 01 xx yy − −∆∆ 03 0 2 1 2 xx yy − −∆∆ 04 0 3 1 3 xx yy − −∆∆ 1 1x 1y 12 12 xx yy − − 13 12 xx yy − −∆∆ 14 1 2 2 2 xx yy − −∆∆ ----------------------- ----------------------- --------------------2 2x 2y 23 23 xx yy − − 24 23 xx yy − −∆∆ ------------------------------------------ ------------------ ----------------------- ----------------------- -------------------- 3 3x 3y 34 34 xx yy − − ------------------------------------------ ------------------ --------------------- --------------------- ------------------ ----------------------- ----------------------- -------------------- 4 4x 4y --------------- -------------------- -------------------- ---------------------- Observação: Para h constante , a relação entre diferenças divididas e finitas é dada por : i k ki k y h!.k 1y ∆∆ = . Exemplo: x f(x) iy∆ i2 y∆ i3 y∆ i4 y∆ 0 0 4 5 1 0 2 8 19 10 1 ----- 3 27 49 14 ----- ----- 5 125 91 ----- ----- ----- 6 216 ----- ----- ----- ----- 39 4-Interpolação Polinomial Exemplo2: A velocidade do som na água varia com a temperatura conforme tabela: 532,1 110 538,1 4,104 544,1 9,98 548,1 3,93 552,1 86 )s/m(velocidade ) C(atemperatur o Solução: 1) Cálculo das diferenças divididas x y 0y∆ 02y∆ 03y∆ 04y∆ 86 1.552 -0.00054794 -0.00001289 -0.114 x 10-5 0.136 x 10-6 93.3 1.548 -0.00071428 -0.00003393 0.213 x 10-5 --------------- 98.9 1.544 -0.00109090 0.175 x 10-5 --------------- --------------- 104.4 1.538 -0.00107142 --------------- --------------- --------------- 110 1.532 --------------- --------------- --------------- --------------- 2) Construção do polinômio: )9.98x)(3.93x)(86x()00001289.0)(3.93x)(86x()00054794.0)(86x(552.1)x(p −−−+−−−+−−+= (-0.114 x 10 5− )+ )4.104x)(9.98x)(3.93x)(86x( −−−− (0.136 x 10 6− ) Simplificando a expressão , encontramos 5036369194.0x830077947032.0x99260000534327.0x10840.13666902p(x) 234-6 −+−×= 3) Verificação da validade de p(x) calculado no item 2: 1.55199999)86(p = ; 548.1)3.93(p = ; 544.1)9.98(p = ; 537999999.1)4.104(p = ; 532.1)110(p = 4) Estimativa da velocidade do som quando a temperatura da água atinge 100 C0 , é 1.54293925)100(p ≅ m/s 6) Análise gráfica do problema: 40 E.B.Hauser – Cálculo Numérico 4.5 Erro de Truncamento [ ]nxxnfn xxE ,0),( )1( )!1( )()( ∈++= ξξ ϕ com i 1n n101n 0 1n n10 y)*xx)...(xx)(xx( h)!1n( f)xx)...(xx)(xx()x( ++ + −−−=+−−−= ∆ ∆ϕ , pois i k ki k f h!.k 1y ∆∆ = , para h constante. 4.6 Aplicações utilizando o sistema Maple APLICAÇAO 1 Cinquenta animais de uma espécie ameaçada de extinção foram colocados numa reserva e um controle populacional mostrou os dados: 69 13 73 10 76 7 77 5 73 3 60 1 50 0 animaisdequantidade )anos(t a) Determinar a matriz de Vandermonde para o problema e determinar o valor do respectivo Número de Condicionamento (Cond e Determinante Normalizado). O que podemos concluir? b) Determinar o polinômio interpolador utilizando todos os dados tabelados. c) Verificar a validade do modelo encontrado. d) Estimar a a população no quarto ano. É possível estimar a população no décimo quinto ano utilizando o polinômio determinado no ítem b. e) Em que ano essa população animal atingiu seu máximo? Qual a população máxima atingida? f) Plotar num mesmo sistemas de eixos os pontos tabelados e a e o polinômio interpolador determinado no ítem b. Resposta: ⎥⎥ ⎥⎥ ⎥⎥ ⎥⎥ ⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎢⎢ ⎢ ⎣ ⎡ = 65432 65432 65432 65432 65432 1313131313131 1010101010101 7777771 5555551 3333331 1111111 0000001 V Determinante Normalizado= 0.9072675023 x 10 -11 Cond (V)= 39124291.11 p(t) = .5095984263e-4*t^6-.2488977072e-2*t^5+.4187474562e-1*t^4-.2409954552*t^3- .6536635972*t^2+10.85522232*t+50. p(4) ≅ 75.91851648 população máxima ≅ p(5.312779138)= 77.05456458 41 4-Interpolação Polinomial APLICAÇÃO 2 Para determinar a resistência elétrica de um solo num sistema de aterramento, enterra-se duas hastes de cobre e aplica-se uma determinada voltagem, resultando numa corrente elétrica. Numa experiência deste tipo, foram obtidos os seguintes dados: 5,4 50 3,4 47 5,3 40 8,2 35 2 30 ))A(amperecorrente(y ))V(voltsvoltagem(x − − Estimar a corrente se a voltagem aplicada for de 43A usando o Polinômio Interpolador de Newton. Resposta: p(43)≅ 3,88 APLICAÇÃO 3 “Ao considerar que no Japão a vida média já é superior a 81 anos, a esperança de vida no Brasil de pouco mais que 71 anos ainda é relativamente baixa. E, de acordo com a projeção mais recente da mortalidade, somente por volta de 2040 o Brasil estaria alcançando o patamar de 80 anos de esperança de vida ao nascer. (Ver www.ibge.gov.br em População / Tábuas Completas de Mortalidade). A esperança de vida ao nascer de 71,3 anos coloca o Brasil na 86ª posição no ranking da ONU, considerando as estimativas para 192 países ou áreas no período 2000-2005 (World Population Prospects: The 2002 Revision; 2003). Entre 1980 e 2003 a esperança de vida ao nascer, no Brasil, elevou-se em 8,8 anos: mais 7,9 anos para os homens e mais 9,5 anos para as mulheres. Em 1980, uma pessoa que completasse 60 anos de idade teria, em média, mais 16,4 anos de vida, perfazendo 76,4 anos. Vinte e três anos mais tarde, um indivíduo na mesma situação alcançaria, em média, os 80,6 anos. Aos 60 anos de idade os diferenciais por sexo já não são tão elevados comparativamente ao momento do nascimento: em 2003, ao completar tal idade, um homem ainda viveria mais 19,1 anos, enquanto uma mulher teria pela frente mais 22,1 anos de vida”.Na tabela acima obtemos informações sobre a esperança de vida às idades exatas, especialmente, a esperança de vida ao nascer que expressa o número médio de anos que se espera viver um recém-nascido (que, ao longo da vida, estivesse exposto aos riscos de morte da tábua de mortalidade em questão 42 E.B.Hauser – Cálculo Numérico http://www.ibge.gov.br/home/presidencia/noticias/noticia_visualiza.php?id_noticia=266&id_pagina=1). Além dos múltiplos usos, não somente pela Demografia e Previdência Provada, mas também pelas demais Ciências Sociais, a tabela é utilizada para determinar, juntamente com outros parâmetros, o chamado fator previdenciário para o cálculo das aposentadorias das pessoas regidas pelo Regime Geral de Previdência Social. TAREFA: Considerar os resultados de 2003. 0.15 70 4.18 65 1.22 60 0.26 55 1.30 50 3.48 30 0.53 25 8.57 20 6.62 15 5.67 10 2.75 0 Mulheres)2003em,anos(vidadeectativaexp 2003emidade − 1.13 70 9.15 65 1.19 60 5.22 55 2.26 50 5.42 30 8.46 25 0.51 20 5.55 15 4.60 10 6.67 0 Homens)2003em,anos(vidadeectativaexp 2003emidade − A) Determinar o polinômio interpolador utilizando todos os dados tabelados. Sugestão: ?interp B) Verificar a validade do modelo encontrado. C) Plotar num mesmo sistema de eixos os pontos tabelados e a e o polinômio interpolador determinado no item b. D)Estimar a expectativa de vida para pessoas com idade em 2003, variando de 14 a 30 anos. > Estimativas_mulheres:= array( [ seq( [i, pexpvidam(i)], i = 14..30)]); > Estimativas_homens:= array( [ seq( [i, pexpvidah(i)], i = 14..30)]); := Estimativas_mulheres ⎡ ⎣ ⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢ ⎤ ⎦ ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥ 14 63.56097596 15 62.6000001 16 61.641461117 60.6831095 18 59.7236187 19 58.7625088 20 57.8000000 21 56.8368317 22 55.8740682 23 54.9129142 24 53.954550 25 53.000000 26 52.050032 27 51.105097 28 50.165308 29 49.230445 30 48.300001 := Estimativas_homens ⎡ ⎣ ⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢ ⎤ ⎦ ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥ 14 56.46550812 15 55.5000000 16 54.5559710 17 53.6353083 18 52.7375232 19 51.8602890 20 51.0000001 21 50.1523108 22 49.3126236 23 48.476507 24 47.640028 25 46.800001 26 45.954138 27 45.101122 28 44.240601 29 43.373125 30 42.500001 43 4-Interpolação Polinomial 4.7 Exercícios (Fonte: Cálculo Numérico Computacional. Claudio, Dalcídio Moraes e Marins, Jussara M. São Paulo: Ed.Atlas,1994.) 1. A tabela abaixo dá o volume de água num tanque elástico (usado para transpor-te de óleo, leite, etc. em caminhões) para várias cotas de água. Determinar y(0,12). x(m) 0,1 0,6 1,1 1,6 2,1 y )3(m 1,1052 1,8221 3,0042 4,9530 8,1662 2. Uma hidroelétrica tem capacidade máxima de 60Mw, a qual é determinada por três geradores de respectivamente 30Mw, 15Mw e 15Mw. A demanda de energia varia num ciclo de 24h e é uma função dela que o engenheiro ope- racional distribui as tarefas dos geradores. Sabe-se que a demanda mínima ocorre entre 1 e 5h e a demanda máxima entre 13 e 17h . Pede-se para achar a partir dos dados abaixo essas demandas máximas e mínimas . H 2 3 4 5 13 14 15 16 17 Demanda (Mw) 16,4 15,2 14,9 16 28 36,5 43 34 31,2 3. Um paraquedista realizou seis salto; saltando de alturas distintas em cada salto, foi testada a precisão de seus salto em relação a um alvo de raio de 5m, de acordo com a altura. A distância apresentada na tabela abaixo é relativa à circunferência. ALTURA (m) DISTÂNCIA DO ALVO 1O SALTO (1500) 2O SALTO (1250) 3O SALTO (1000) 4O SALTO (750) 5O SALTO (500) 35 25 15 10 7 Levando em consideração os dados acima, a que provável distância do alvo cairia o paraquedista se ele saltasse de uma altura de 850m ? 4. Um veículo de fabricação nacional, após vários testes, apresentou os resultados abaixo, quando se analisou o consumo de combustível de acordo com a velocidade média imposta ao veículo. Os testes foram realizados em rodovia em operação normal de tráfego, numa distância de 72 km. Verificar o consumo aproximado para o caso de ser desenvolvida a velocidade de 80km/h. VELOCIDADE (km/h) CONSUMO (km/l) 55 14,08 70 13,56 85 13,28 100 12,27 120 11,3 140 10,4 44 E.B.Hauser – Cálculo Numérico 5. Numa esfera de superfície conhecida, o coeficiente de absorção 0,7 foi mantido à temperatura de 6000 Ko . Foi calculada a energia irradiada de acordo com o tempo de irradiação, obedecendo à tabela . Pede-se para obter a possível energia irradiada quando a irradiação atingir o tempo de 25 minutos. ENERGIA IRRADIADA (Joules) TEMPO DE IRRADIAÇÃO (s) 71,72 . 310 600 94,72 . 310 800 118,4 . 310 1000 142,08 . 310 1200 165,76 . 310 1400 189,44 . 310 1600 6. Uma corda foi tensionada sob a ação de pesos distintos, quando para os respectivos pesos foram calculadas as devidas velocidades de propagação que estão indicadas abaixo. Pede-se para calcular a velocidade de propagação quando a corda está tensionada sob a ação de um peso de 7250 gf. PESO (gf) VELOCIDADE (cm/s) 6000 13728,13 6500 14288,69 7000 14828,07 7500 15348,51 8000 15851,87 Respostas 1) 126904937,1)12,0(p ≅ 2) A demanda mínima é 8632,14 Mw. e ocorre entre 3h e 4h da manhã. A demanda máxima é 101,43 Mw. e ocorre entre 14h e 15h. 3) 4128,11)850(p = 4) 4512685,13)80(p 46701783,13)80(p = = 5) 5937,177618)1500(p = 6) p(7250) = 15090,53234 45 5. Ajuste de Funções Aplicação1: Os dados acima tabelados descrevem a intensidade da luz como uma função da distância da fonte, I(d), medida num experimento. Determinar I (d) ≅ Y (d) = CBdAd 1 2 ++ . 15.0 75 18.0 70 21.0 65 24.0 60 28.0 55 34.0 50 42.0 45 52.0 40 67.0 35 85.0 30 I d Aplicação 2:Segundo a lei de Boyle, o volume de um gás é inversamente proporcional à pressão exercida, mantendo-se constante a temperatura. Para um certo gás, foram observados os seguintes valores: 04,1 5,3 15.1 0,3 28,1 5,2 43,1 0,2 60,1 5,1 91,1 0,1 27,2 5,0 85,2 0,0 Volume essãoPr Ajustar os dados tabelados a uma hipérbole do tipo: Vol(p) ≅ pB A)p(Y += 46 E.B.Hauser – Cálculo Numérico 5.1 Introdução O ajustamento é uma técnica de aproximação. Conhecendo-se dados experimentais , , )y,x( ii , n,...,1,0i = , deseja-se obter a lei )x(fy = relacionando x com y. Devido aos erros experimentais nos n+1 pares, ))x(f,x( ii , teremos em geral 11 )x(f ε+ ; 22 )x(f ε+ ; ... ; nn )x(f ε+ , isto é , é impossível calcular exatamente a função f(x). Por isso , em vez de procurarmos a função f tal que passa por cada um dos pontos experimentais, determinaremos a função que melhor se ajusta aos pontos dados. O ajustamento traduz um comportamento médio. Para ajustar uma tabela de dados a uma função devemos: • Conhecer a natureza física do problema ; • Determinar o tipo de curva a que se ajustam os valores tabulados (graficamente e/ou cálculo das diferenças finitas ou divididas) ; • Calcular os parâmetros da curva. 5.2. Escolha da Função de Ajuste a) Função Linear (regressão linear) : xaa)x(Y 10 += , se iy∆ ≅ cte ou iy∆ ≅ cte . b) Parábola (ajuste quadrático): 2210 xaxaa)x(Y ++= , se iy2∆ ≅ cte ou .2 cteiy ≅∆ c) Polinômio de grau p: se i p y∆ ≅ cte ou .cteyip ≅∆ d) Função Exponencial: xab)x(Y = , se .cte x ylog i i ≅∆ ∆ d) FunçãoPotência: pax)x(Y = , cte xlog ylog i i ≅∆ ∆ Outros tipos de Funções Ajuste: • xaa y 1 xaa 1)x(Y 10 10 +=⇒+= ; ( ) ( ) .cte x y/1y/1 i i i ≅= ∆ ∆∆ • ;xaa y x xaa x)x(Y 10 10 +=⇒+= ( ) ( ) ( ) .ctex y/xy/x i ii ii ≅= ∆ ∆∆ • 47 5 - Ajuste de Funções • 22102 210 xaxaa Y 1 xaxaa 1)x(Y ++=⇒++= ; ( ) ( ) ( ) .cte xx y/1y/1 y/1 i2i i1i i 2 ≅− −= + + ∆∆∆ • 2cxbxae)x(Y += ; 2cxbxalnyln ++= , .cte xx ylnyln i2i i1i ≅− − + + ∆∆ Exemplo: Considerando a tabela abaixo, como i,1yi 3 ∀=∆ , é adequado o ajuste dos dados abaixo tabelados a um polinômio de grau 3. i xi f(xi)=yi iy∆ iy2∆ iy3∆ iy4∆ 0 0 0 4 5 1 0 1 2 8 19 10 1 ----- 2 3 27 49 14 ----- ----- 3 5 125 91 ----- ----- ----- 4 6 216 ----- ----- ----- ----- 5.3 - Determinação dos Parâmetros da Função de Ajuste CRITÉRIO DOS MINÍMOS QUADRADOS Seja pxpaxaxaaY ++++= ...2210 a função de ajustamento. Dada uma tabela com n+1 pontos ( )iyix , , chamamos resíduo a diferença entre o valor de iY da equação de ajustamento e o valor tabulado de iy . n,...,1,0i,yY iii ==− δ . O critério dos mínimos quadrados estabelece que: ∑ = = n i mínimoi 0 2δ . Seja ()∑ = −= n i iyiYpaaaF 0 2),...1,0( . Para F ter valor mínimo , é preciso que 0 a F 0 =∂ ∂ ; 0 a F 1 =∂ ∂ ; . . . ; 0=∂ ∂ pa F 48 E.B.Hauser – Cálculo Numérico 5.31 -Ajuste Linear A função de ajuste terá a forma xaa)x(Y 10 += . Pelo critério dos Mínimos Quadrados é necessário que : ∑ ∑ = = −+=−= n i n i iyixaaiyiYF 0 0 2)10( 2)( deve ser mínimo . Sendo F uma função de duas variáveis, 0a e 1a , o menor valor de F será obtido através de : 0 0 = a F δ δ ; 0 1 = a F δ δ e assim: 0 0 = a F δ δ , ( )∑ = =−+ n i iyixaa 0 0102 e 0 1 = a F δ δ , ( )∑ = =−+ n i ixiyixaa 0 0102 . Construímos o sistema de duas equações e duas variáveis: ( ) ( ) ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ = =−+ = =−+ ∑ ∑ n i ixiyxaa n i iyixaa 0 01102 0 0102 ou ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ = = = += = = = += ∑ ∑ ∑ ∑ ∑ ∑ n i n i n i ixaixaiyix n i n i n i ixaaiy 0 0 0 2 10 0 0 0 10 . Obtemos: ( ) ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ = = = =+ = = =++ ∑ ∑ ∑ ∑ ∑ n i n i n i iyixixaixa n i n i iyixaan 0 0 0 2 10 0 0 101 Resolvendo-se este último sistema linear são obtidos os valores de 0a e 1a e assim determina-se a função de ajuste : xaaY 10 += . 49 5 - Ajuste de Funções Exemplo2 : Dada a tabela , achar a equação da reta que se ajusta usando o método dos Mínimos Quadrados. i ix iy iyix 2ix i Y ( )2iyiY − 0 0 2 0 0 2.07142857 0.005102041 1 1 3 3 1 3.14285714 0.02040816 2 2 5 10 4 4.21428571 0.61734694 3 3 5 15 9 5.2857143 0.08163266 4 4 5.5 22 16 6.3571429 0.73469755 5 5 8 40 25 7.4285714 0.32653061 ∑ 15 28.5 90 55 --------------- 1.78571440 Seja xaaxY 10)( += a função que ajusta os dados . Os parâmetros 0a e 1a constituem a solução do sistema : ∑ ∑ ∑ ∑ ∑ += ++= 2 10 10)1( ixaixaiyix ixaaniy ⇒ ⎩⎨ ⎧ =+ =+ 90155015 5.2811506 aa aa ⇒ ⎩⎨ ⎧ = = 071428571.11 071428572.20 a a Logo, xY 071428571.1071428571.2 += . 50 E.B.Hauser – Cálculo Numérico 5.3.2 - Ajuste Quadrático Seja 2210 xaxaaY ++= a função de ajuste. Pelo critério dos Mínimos Quadrados : ( ) ( )∑ ∑ −++=−= = = n 0i n 0i 2 i 2 i2i10 2 ii210 yxaxaayY)a,a,a(F assume valor mínimo se 0 a F a F a F 210 === δ δ δ δ δ δ . Assim, os parâmetros 210 aea,a são obtidos resolvendo: ( ) ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ ++= ++= +++= ∑ ∑ ∑∑ ∑ ∑ ∑ ∑ ∑ ∑∑ 4 ix2a 3 ix1a 2 ix0aiy 2 ix 3 ix2a 2 ix1aix0aiyix 2 ix2aix1a0a1niy Exemplo3 : Encontre a expressão do polinômio de 2o grau que se ajusta aos dados da tabela abaixo. i ix iy iyix 2ix iyix 2 3ix 4 ix i y∆ iy 2∆ 0 -2 -0.01 0.02 4 -0.04 -8 16 0.52 -0.21 1 -1 0.51 -0.51 1 0.51 -1 1 0.31 -0.25 2 0 0.82 0 0 0 0 0 0.06 -0.13 3 1 0.88 0.88 1 0.88 1 1 -0.07 -0.25 4 2 0.81 1.62 4 3.24 8 16 -0.32 ---------- 5 3 0.49 1.47 9 4.41 27 81 ---------- ---------- ∑ 3 3.5 3.48 19 9 27 115 ---------- ---------- ⎪⎩ ⎪⎨ ⎧ =++ =++ =++ 92115127019 48.322711903 5.32191306 aaa aaa aaa ⇒ 806285714.0201.02102142857.0)( ++−= xxxY 51 5 - Ajuste de Funções 5.3.3-Ajuste a polinômio de Grau p O ajuste a um polinômio de grau p. nppxpaxaxaaY <++++= ,...2210 , exige resolver o sistema linear: ⎥⎥ ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎣ ⎡ ∑ ∑ ∑ ∑ ∑ ∑ ∑ ∑ ∑ ∑ ∑∑+ + + p2 i 1p i p i 1p i 4 i 3 i 2 ii p i 3 i 2 ii x.................xx .................................................. x...xxxx x...xxx)1n( ⎥⎥ ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎣ ⎡ pa a a ..... ..... 1 0 = ⎥⎥ ⎥⎥ ⎥⎥ ⎦ ⎤ ⎢⎢ ⎢⎢ ⎢⎢ ⎣ ⎡ ∑ ∑ ∑ iy p ix iyix iy ........ Exemplo5 – Utilizando o sistema Maple > xv:=[0,2,3,5,6,8]; > yv:=[p(0),p(2),p(3),p(5),p(6),p(8)]; > g:=fit[leastsquare[[x,y],y=a*x^3+b*x^2+c*x+d,{a,b,c,d}]]([xv,y v]); > gll :=unapply(rhs(g),x): > plots[display]({plots[pointplot]([0,p(0),2,p(2),3,p(3),5,p(5), 6,p(6),8,p(8)]),plot(gll(x), x=-.5..8.5, y=- 150..280),plot(gll(x), x=-.5..8.5, y=-150..280,thickness=2)}); 52 E.B.Hauser – Cálculo Numérico 5.3.4 -AJUSTE NÃO LINEAR NOS PARÂMETROS: CASOS REDUTÍVEIS AO LINEAR OU PARABÓLICO POR MUDANÇA DE VARIÁVEIS 5.3.4.1 – Ajuste por Função Exonencial Seja xabY = a função de ajuste. Linearizamos Y, aplicando log ou ln (escolher a base adequada) { { { 1xa0ay blnxalnYln += Determinamos a e b resolvendo: ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ = = = + = = = = ++ ∑∑∑ ∑∑ iyln n 0i ix n 0i 2 ixbln n 0i ixaln n 0i iyln n 0i ixblnaln)1n( Então ⎪⎩ ⎪⎨ ⎧ =→= =→= 1a1 0a0 ebblna eaalna Exemplo 4: Ajustar os dados abaixo a uma função exponencial do tipo xabY = . i ix iy iyln iyix ln 2ix i Y 2)( iyiY − 0 0 3 1.09861 0 0 2.98422 0.0002490 1 0.5 4 1.38629 0.693145 0.25 4.2228031 0.0496412 2 1 6 1.79176 1.79176 1 5.9754291 0.0006025 3 1.5 9 2.19722 3.29583 2.25 8.4555298 0.2964477 4 2 12 2.48491 4.96982 4 11.964948 0.0012286 5 2.5 17 2.83321 7.08303 6.25 16.930930 0.0047706 6 3 24 3.17805 9.53415 9 23.958013 0.0017628 7 3.5 33 3.49651 12.2378 12.25 33.901647 0.8129689 8 4 48 3.87120 15.4848 16 47.972329 0.0007656 ∑ 18 ------------- 22.3378 55.0903 51 ------------- 1.1684403 Os gráficos a seguir ilustram o efeito da linearização dos dados. 53 5 - Ajuste de Funções ⎩⎨ ⎧ =+ =+ 0903.55151018 3378.2211809 aa aa ⇒ 984218125.20093337778.11 002347015.2169432.00 ==→= ==→= aeaa aeba ⇒ xY )00235.2(98422.2= - 5.3.4.2 -AJUSTE POR FUNÇÃO POTÊNCIA Seja baxY = a função de ajuste. Linearizando a função Y(x) , temos : .lnlnln xbay += Resolve-se o sistema de equações lineares e encontra-se a e b: 1 1 0 1 xxln ab aaln yyln = = = = ⇒ ⎪⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ = = = + = = = = ++ ∑∑∑ ∑∑ n 0i iylnixln n 0i 2)ix(lnb n 0i ixlnaln n 0i iyln n 0i ixlnbaln)1n( então 10 a 0 abeeaaaln ==⇒= . Exercício : Os dados abaixo dão a duração de uma broca em função da velocidade de corte. Pede-se para fazer uma tabela de D=D(v) para 100 (10) 180. 8.2 180 9.7 150 28 120 79 100 .)( )/( segD smV D(v)=0.813947103e14*1/(x^5.680591334) - 54 E.B.Hauser – Cálculo Numérico Efeito da linearização dos dadoss 5.3.4.3 -AJUSTE POR FUNÇÃO HIPERBÓLICA : xaa 1Y 10 + = . Linearizando, temos : ⎪⎪⎩ ⎪⎪⎨ ⎧ ∑=∑+∑ ∑ ∑=++ ⇒+= === = = n 0i ii
Compartilhar