Mudanças entre as edições de "Map0151"

Edição atual tal como às 22h23min de 4 de janeiro de 2014

[editar] Números no Computador

Muitos dos algoritmos de cálculo numérico, por uma questão prática, deverão ser executados numa máquina real. Nestas máquinas, no entanto, a capacidade de memória para representação dos números é finita. Vários números reais (infinitos, de fato) terão a mesma representação no computador (ou calculadora), daí originando-se os erros de arredondamentos. Vamos ver qual é a técnica usada atualmente para diminuir os erros de arredondamentos.

[editar] Representação de números inteiros numa base

Seja $\beta > 1$ um número natural. Então podemos representar qualquer número inteiro $k$ como a soma $$k= \text{sgn}(k)*a_0\beta^0 + \cdots + a_s \beta^s$$ onde cada algarismo $a_j$ representa um número natural entre $0$ e $\beta -1$. Só $a_s$ é diferente de zero. Esta representação é única. A representação de $k$ na base $\beta$ é $$k=\pm [a_s\dots a_0]_{\beta}$$ Deixando o sinal de $k$ de fora, para simplificar. Podemos calcular os algarismos do número daddo na base $\beta$ aplicando repetidas vezes o algoritmo da divisão. Escrevendo $$k = a_0 + \beta(a_1 +\beta( a_2 + \cdots +\beta a_s)\cdots ))$$ vemos que $a_0$ é o resto da divisão $k//\beta$, $a_1$ é o resto da divisão $(k-a_0)//\beta$ e assim por diante.

Exemplo: $39 = [100111]_2$

[editar] Representação de números fracionários e decimais numa base

Seja novamente $\beta > 1$ um número natural.

Se $x \in (0,1)$ então $$x= \frac{b_1}{\beta} + \cdots + \frac{b_k}{\beta^k} + \cdots$$ Diremos que $x=[0.b_1 b_2\dots]_{\beta}$ é a representação fracionária na base $\beta$ de $x$. Quando não houver dúvidas de que base se trata, omite-se a base da notação!

Exemplo: $0.9= [0.1110011001100...]_2$

[editar] Uma função para colocar um número decimal na forma binária

[editar] Representação em ponto flutuante

Consideramos uma base fixa $\beta > 1$. Um número real $\alpha \in \mathbb{R}$ positivo pode ser escrito nesta base como: $$\alpha = [a_k\cdots a_0.b_1b_2\cdots]_{\beta}$$ Isto significa que: $$\alpha = a_k\beta^k + \cdots + a_0 + \frac{b_1}{\beta} +\frac{b_2}{\beta^2} + \cdots$$ Na equação acima, colocando $\beta^{k+1}$ em evidência temos: $$\alpha = \left( \frac{a_k}{\beta} + \cdots + \frac{a_0}{\beta^{k+1}} + \frac{b_1}{\beta^{k+2}} +\frac{b_2}{\beta^{k+3}} + \cdots\right)\times \beta^{k+1}$$ ou ainda, lembrando da notação de um número numa base dada: $$\alpha = [0.a_k \cdots a_0 b_1 b_2\cdots]_{\beta} \times \beta^{k+1}$$ Esta última fórmula é importante. O número real $\alpha$ fica caracterizado por três dados:

• O número $m = [0.a_k \cdots a_0 b_1 b_2\cdots]_{\beta} \in [0.1,1)$ chamado de mantissa.
• O número $e = k+1$ chamado de expoente
• O sinal do número $\sigma$

Esta representação do número $\alpha$ como $\sigma m \times \beta^e$ chamaremos de representação normal em ponto flutuante na base $\beta$. Em geral a base fica clara pelo contexto!

[editar] Números de máquina

Continuamos com a base fixa ($\beta$), mas na representação normal em ponto flutuante vamos admitir apenas números com a mantissa limitada a $D$ dígitos e o expoente limitado, em módulo, por um número inteiro $M$. O conjunto $$\mathbf{M} = \left\{ [0.d_1\cdots d_D]_{\beta}\times \beta^e: |e| \leq M \right\}$$ é um conjunto de números de máquina definido pelos números $(\beta, D, M)$. Este é um conjunto finito com $(\beta-1)(\beta^{(D-1)})(2M+1)$ elementos. Vamos agora definir a função arredondamento, que é uma função sobrejetora $\operatorname{rd}: \mathbb{R} \to \mathbf{M}$. Seja $\alpha$ um número real; se o expoente de $\alpha$ for maior que $M$ não representamos o número. Senão tomamos sua mantissa e se o dígito $b_{d+1}\geq \beta/2$ escolhemos o menor número de máquina maior que $\alpha$ para representá-lo. No caso de $b_{d+1} <\beta/2$ o representante será o maior número de máquina menor que $\alpha$.

Se denotarmos por $\underline{\alpha}$ o truncamento do número $\alpha$, isto é, o maior número de máquina menor que $\alpha$. E por $\overline{\alpha}$ o menor número de máquina maior que $\alpha$. Temos que: $$\overline{\alpha} - \underline{\alpha} = \beta^{e-D}$$ Com isto temos uma avaliação do erro absoluto e erro relativo do arredondamento: $$|\alpha -\operatorname{rd}(\alpha)| \leq \frac{\beta^{e-D}}{2}$$ $$\frac{|\alpha -\operatorname{rd}(\alpha)|}{|\alpha|}\leq \frac{\beta^{1-D}}{2}$$

[editar] Aritmética de números em pontos flutuantes

A extensão das operações aritméticas de soma e multiplicação de números reais faz-se pelo arredondamento: $$\alpha \oplus \beta = \operatorname{rd}( \operatorname{rd}(\alpha) + \operatorname{rd}(\beta))$$ $$\alpha \otimes \beta = \operatorname{rd}( \operatorname{rd}(\alpha) \times \operatorname{rd}(\beta))$$ Estas operações deixam de ser associativas.

[editar] Aproximações usando polinômios de Taylor

Para fazer as avaliações de funções de variável real, utilizaremos apenas as operações aritméticas. Basicamente, isto significa que só usaremos polinômios. Em primeiro lugar veremos o método de avaliação de Horner.

[editar] Método de Horner

É a forma de se avaliar o polinômio num ponto dado executando o menor número de operações aritméticas. O polinômio da forma $$p(x) =a_nx^n+ a_{n-1}x^{n-1}+\cdots + a_0$$ no ponto $x_0$ é avaliado na ordem $$p(x)= [[[[a_nx_0 + a_{n-1}]x_0+a_{n-2}]x_0+\cdots]x_0 +a_0]$$ Isto é, calculamos a sequência

• $A_0 = a_n$
• $A_{k+1} = A_kx_0+ a_{n-k}$ até o $A_n$ que é a avaliação.

[editar] Teorema de Taylor

Teorema: Seja $f:(a,b) \to \mathbb{R}$ uma função de classe $\mathcal{C}^r$, talque exista a derivada $f^{(r+1)}(c), \forall c \in (a,b)$. Então para todo $x$ do intervalo, e $c$ fixado existe $\xi \in |c,x|$ tal que: $$f(x) = f(c) + f^\prime(c)(x-c)+ \cdots + \frac{f^{(r)}(c)}{r!}(x-c)^{r} + \frac{f^{(r+1)}(\xi)}{(r+1)!}(x-c)^{r+1}$$

[editar] Zeros de funções

Dada uma função $f:(a,b)\to \mathbb{R}$, queremos encontrar um $r\in (a,b)$ tal que $f(r)=0$. Um tal número será chamado de zero da função $f$. Para assegurar a existência de zeros no intervalo, precisamos de alguma regularidade da função. Se $f$ é uma função contínua, uma condição suficiente para a existência do zero é que $f(a)*f(b) < 0$, e uma técnica para chegar ao zero é ir contraindo este intervalo, mantento a relação $f(a)*f(b) < 0$, até chegarmos suficientemente próximo. No método da bissecção vamos dividindo o intervalo ao meio. Se $(a_n,b_n)$ é o intervalo com $f(a_n)*f(b_n) < 0$, definimos $x_{n+1} = (a_n+b_n)/2$ e o novo intervalo como:

• $(a_{n+1},b_{n+1})=(x_{n+1},b_n)$ se $f(a_n)*f(x_{n+1}) < 0$
• $(a_{n+1},b_{n+1})=(a_n,x_{n+1})$ caso contrário

O tamanho do intervalo no n-ésimo passo é $|b-a|/2^n$. A seguir uma função em python que implementa a bissecção.

[editar] Algoritmo da bissecção

 
|  k  | x0       | x1       | xmid     |f(x0)*f(xmid)|
| 0   | 0.000000 | 2.000000 | 1.000000 | +2.00000000 |
| 1   | 1.000000 | 2.000000 | 1.500000 | -0.25000000 |
| 2   | 1.000000 | 1.500000 | 1.250000 | +0.43750000 |
| 3   | 1.250000 | 1.500000 | 1.375000 | +0.04785156 |
| 4   | 1.375000 | 1.500000 | 1.437500 | -0.00726318 |
| 5   | 1.375000 | 1.437500 | 1.406250 | +0.00245667 |
| 6   | 1.406250 | 1.437500 | 1.421875 | -0.00048804 |
| 7   | 1.406250 | 1.421875 | 1.414062 | +0.00000960 |
| 8   | 1.414062 | 1.421875 | 1.417969 | -0.00000454 |
| 9   | 1.414062 | 1.417969 | 1.416016 | -0.00000218 |
| 10  | 1.414062 | 1.416016 | 1.415039 | -0.00000100 |
| 11  | 1.414062 | 1.415039 | 1.414551 | -0.00000041 |
| 12  | 1.414062 | 1.414551 | 1.414307 | -0.00000011 |
| 13  | 1.414062 | 1.414307 | 1.414185 | +0.00000004 |

[editar] Algoritmo de Newton

O algoritmo de Newton-Raphson usa propriedades de diferenciabilidade da função $f$ para criar uma sequência recorrente que converge ao zero da função que atrairá um ponto $x_0$, que é uma estimativa inicial. A fórmula para gerar a sequêmcia é: $$x_{k+1} = x_k - \frac{f(x_k)}{f^\prime(x_k)}$$ Abaixo o algoritmo em python:

| k | x_k       | |x_k-xk-1|   |
|1  | 1.4166667 | 0.0833333
|2  | 1.4142157 | 0.0024510
|3  | 1.4142136 | 0.0000021

[editar] Algoritmo da secante

O método da secante pode ser visto como uma adaptação do método de Newton, onde a derivada da função $f$ é computada pela fórmula de aproximação $$f^\prime(x_k) \sim \frac{f(x_k)-f(x_{k-1})}{(x_k - x_{k-1})}$$ definindo a sequência $$x_{k+1} = x_k - \frac{f(x_k)(x_k -x_{k-1})}{(f(x_k)-f(x_{k-1}))}$$

O código em Python pode ser o seguinte

| k | x_k       | f(xk)   |
|1  | 2.0000000 | 2.0000000
|1  | 1.0000000 | -1.0000000
|2  | 1.3333333 | -0.2222222
|3  | 1.4285714 | 0.0408163
|4  | 1.4137931 | -0.0011891
|5  | 1.4142114 | -0.0000060

[editar] Resolução de sistemas lineares

Do ponto de vista de algoritmos numéricos existem duas formas para se resolver um sistema de equações lineares da forma: $$\begin{align} a_{11}x_1 + \cdots + a_{1n}x_n &=& b_1 \\ \vdots & & \\ a_{n1}x_1 + \cdots +a_{nn}x_n &= &b_n \end{align}$$

Os Métodos Diretos têm, em geral, duas fases. Numa primeira fase, um sistema genérico é transformado, numa sequência finita de passos, num sistema de uma classe de sistemas, do qual se conhece um algoritmo eficiente para a solução. Na segunda fase, resolvemos o sistema gerado. O método da eliminação de Gauss (MEG) é um método direto que descreveremos a seguir.

Os Métodos Iterativos são obtidos por uma sequência de estimativas da solução. Esta sequência deve convergir para a solução, e muito raramente se obtem a solução exata. Mas tem-se um bom controle de erros. Descreveremos o método de Gauss-Seidel que ilustra esta técnica.

[editar] Método da Eliminação de Gauss

Em primeiro lugar vamos descrever os algoritmos para a resolução de sistemas triangulares. Uma matriz quadrada $U = (u_{ij})$ é triangular superior quando tivermos a condição $u_{ij} = 0$ para $i> j$. Neste caso podemos resolver o sistema linear $Ux=b$,usando o algoritmo de baixo pra cima em $n$ passos:

• Passo 1 : Encontro $x_n = b_n/u_{nn}$
• Passo k : Encontro $x_{n-k+1} = \frac{1}{u_{n-k+1 n-k+1}}\left(b_{n-k+1} - \sum_{j=n-k+2}^{n} u_{n-k+1 j}x_j\right)$

A matriz $L=(l_{ij})$ é triangular inferior quando $l_{ij}=0 \text{ e } l_{ii}=1$ para $i < j$. Para resolver os sistema $Lx=b$ usando um algoritmo parecido com o anterior.

• Passo 1 : Encontro $x_1 = b_1/l_{11}$
• Passo k : Encontro $x_{k} = \frac{1}{u_{kk}}\left(b_{k} - \sum_{j=1}^{k-1} u_{k j}x_j\right)$

No método da eliminação de Gauss (MEG), vamos transformar um sistema linear numa sistema equivalente triangular superior, isto é, quando a matriz dos coeficientes $A$ é triangular superior. Faremos isso usando as operações elementares, que descrevemos a seguir:

• $\mathcal{O}(i,j)$ é a troca da linha $i$ com a linha $j$
• $\mathcal{O}(i,\alpha)$ é a multiplicação da linha $i$ por um número, não nulo, $\alpha$
• $\mathcal{O}(i,\alpha, j)$ é a substituição da linha $i$ pela soma desta com $\alpha$ multiplicado pela $i-$ésima linha.

O algoritmo é a aplicação sucessivas dos seguintes passos com $k$ variando de $1$ até $n-1$

 Passo $k$: Caso $a_{kk} \neq 0$, faça de $i=k+1$ até $i=n$ $$m_{ik} = a_{ik}/a_{kk}$$  e aplique ao sistema a operação elementar $\mathcal{O}(i,m_{ik}, k)$.

Ao fim do passo $n-1$ o sistema estará triangularizado. Para evitar a possibilidade de não se triangularizar o sistema, mesmo ele sendo possível e determinado, executamos a pivotação que é o seguinte procedimento: antes da execução do passo $k$ do MEG. Buscamos, na $k-$ésima coluna, a partir da $k-$ésima linha o índice $i$ tal que $|a_{ik}|$ seja o maior possível, e executamos $\mathcal{O}(i,k)$. Com isto garantimos a triangularização de qualquer sistema que tenha solução única.

A seguir o código em python para a eliminação de Gauss

Versão sem pivotação (não sei se precisa!)

[editar] Decomposição LU

Diremos que uma matriz $L$ é triangular inferior se $L=(l_{ij})$ com $$l_{ij}=0 \text{ se } i \lt j \\ l_{ii}=1$$

Por outro lado, a matriz $U$ será triangular superior se $U=(u_{ij})$ com $$u_{ij}=0 \text{ se } j \lt i$$ Uma matriz quadrada $A$ tem uma decomposição LU se existem uma matriz triangular inferior, $L$ e uma triangular superior $U$ tal que temos a igualdade: $$LU=A$$

Uma matriz inversível tem decomposição LU se e somente se seus menores principais têm determinantes não nulos. Neste caso pode-se resolver o sistema $$Ax=b$$ Resolvendo primeiro $$Lz=b$$ usando substituição de cima pra baixo, e depois, resolver $$Ux =z$$ usando substituição de baixo pra cima.

[editar] Métodos Iterativos

A idéia dos métodos iterativos é produzir uma sequência de estimativas vetoriais $\mathbf{x}^{(k)}$, que idealmente convirja para a solução $\mathbf{x}$ e cada passo é calculado a partir da estimativa anterior. Costumam ser bem mais rápidos no caso de matrizes muito grandes. Para construir recursivamente a sequência que dá origem ao método iterativo fazemos uma decomposição da matriz $A$ de coeficientes do sistema como $$A = N - P$$ Onde $N$ é de uma classe de sistemas inversíveis e fáceis de resolver. Então resolver o sistema: $$A\mathbf{x} =b$$ é equivalente a resolver a equação: $$N\mathbf{x} = b+ P\mathbf{x} \text{ ou } \mathbf{x}=N^{-1}(b+ P\mathbf{x})$$ A recursão $$N\mathbf{x}^{(k+1)} = b+ P\mathbf{x}^{(k)}$$ é a que converge para a solução caso $||N^{-1}P|| < 1$ onde a norma de matrizes é induzida de uma norma em $\mathbb{R}^n$. Destacamos dois casos

• Método de Jacobi $$N = \text{diag}(a_{11}, \cdots, a_{nn})$$
• Método de Gauss-Seidel $$N = (n_{ij}) \text{ com } n_{ij} = a_{ij} \text{ se } i\geq j \\ n_{ij}=0 \text{ se } i\lt j$$

[editar] Método dos Mínimos Quadrados: Caso discreto

Consideremos uma tabela de pontos $T = \{(x_0,y_0),\cdots,(x_k,y_k)\}$. Podemos interpretar esta tabela de duas formas: 1. Pode ser um conjunto de informações aproximadas sobre uma função de determinado tipo. E neste caso procuraremos a função daquele tipo que melhor se ajuste à tabela. 2. Pode ser um conjunto de informações precisas sobre uma função diferenciável. Mas não sabemos mais o tipo da função.

No primeiro caso faremos um ajuste da função para a tabela pelo Método dos Mínimos quadrados. A função ajustada não respeita, necessariamente, a tabela dada. Para o segundo caso estudaremos a Interpolação polinomial

Comecemos com o método dos mínimos quadrados (MMQ). Lembremos que temos dadas uma tabela de pontos $T = \{(x_0,y_0),\cdots,(x_k,y_k)\}$ e uma família de funções $\mathcal{F}$ onde o domínio de cada uma das funções de $\mathcal{F}$ contém os pontos $x_i$ da tabela $T$. O resíduo quadrático de uma função $f$ com relação à tabela é: $$r(f,T)= \sum_{i=0}^k \left(y_i -f(x_i)\right)^2$$ O problema de MMQ é encontrar uma função $f^* \in \mathcal{F}$ que minimize $r(f,T)$