AED - Cálculo Número - RelatorioSplineCubico

Prévia do material em texto

PONTIFÍCIA UNIVERSIDADE CATÓLICA DE GOIÁS
FUNDAMENTOS DA COMPUTAÇÃO IV
SPLINE CÚBICA
GOIÂNIA
2015
AED(Atividade Extra Disciplinar) referente 
à 1N2. Relatório sobre Spline Cúbica. Profº 
Gustavo Siqueira Vinhal. Apresentado pela 
aluna Louíse Lima Torrejon, com matrícula 
20131002804960. 
INTRODUÇÃO
A interpolação é um método utilizado para estimar valores de uma
função dentro de um intervalo conhecido, isto é, tendo um conjunto de dados
tabelados, podemos inserir pontos entre os mesmos. Conhecendo os
respectivos conjuntos de dados ou seja, (n+1)
pontos, podemos obter o polinômio interpolador que por sua vez passa
por estes mesmos pontos. O grau de será no máximo n, porém
dependerá do número de pontos que se pretende interpolar e da precisão que
se almeja(RUGGIERO e LOPES, 1997; SPERANDIO, MENDES e SILVA,
2003).
Ruggiero e Lopes evidenciam o papel da interpolação: Interpolar uma
função que consiste em aproximar essa função por outra função 
escolhida entre uma classe de algumas propriedades. A função é então
usada em substituição à função A necessidade de se efetuar esta
substituição surge em várias situações, como quando são conhecidos somente
valores numéricos da função para um conjunto de pontos e é necessário
calcular o valor da função em um ponto não tabelado (RUGGIERO e LOPES,
1997).
A interpolação é aplicável na resolução de problemas como: crescimento
populacional, na resolução de problemas físicos, computacionais e entre
outros.
O grau do polinômio interpolador depende muito do número de pontos
conhecidos e do número de pontos que se quer interpolar. Porém, um
polinômio de grau elevado pode causar extrapolação dos valores de no
intervalo interpolação. Por isso, às vezes por mais que seja o número n de
pontos conhecidos, o polinômio interpolador deve ser obtido usando poucos
pontos conhecidos na vizinhança do ponto que se pretende interpolar.
A figura 1 ilustra uma situação em que o polinômio interpolador extrapola
os valores de 
Um método de interpolação que usa poucos pontos para obter o
polinômio interpolador é a interpolação por Spline, sendo esta função denotada
por onde seu grau é p, tendo os pontos com i = 0, 1, …, n.
Figura 1 – Função e seu polinômio interpolador 
Como observado na Figura 1, nas situações em que o número de pontos
de interpolação e grande, a aproximação obtida com um polinômio de grau
elevado é dominada pelos erros de arredondamento. Sendo que isso nos
sugere uma interpolação por Spline, ou seja, uma interpolação por partes.
FUNÇÃO SPLINE EM INTERPOLAÇÃO
Os primeiros estudos sobre Splines foram feitos por Euler (1707-1783) e
outros autores, mas no início do século XX R. Courant (1943) e I. Shoemberg
(1946) formularam estudos mais detalhados sobre Splines, que chegaram aos
dias de hoje. Em geral, a teoria de splines foi desenvolvida a partir das
necessidades práticas de aproximação (SCHEID, 1991).
Na aproximação de uma função por uma função linear por partes:
A spline cúbica, usa polinômios cúbicos para interpolar os nós dois a
dois.
Características positivas:
* Flexibilidade do polinômio cúbico.
* Sem picos ou trocas de curvatura abruptas nos nós.
Supondo que esteja tabelada em (n+1) pontos, o polinômio
interpolante é composto de n polinômios , um para casa intervalo
entre dois nós. 
As condições de existência dos polinômios são:
i) para com k variando de 1 a n;
ii) com i variando de 0 a n;
iii) k de 0 a (n-1);
iv) 
v) 
Os polinômios k de 1 a n, são trabalhados a partir da seguinte
expressão, de forma a simplificar os cálculos:
* 4n coeficientes;
* 4n equações.
Número de equações obtidas pelas condições de interpolação:
* (n+1) para interpolar nos nós;
* (n-1) para ser contínua nos nós;
* (n-1) para as derivadas de nos nós;
* (n-1) para as derivadas segundas de nos nós.
Total: 4n-2 equações → Insuficientes
Nesta análise as duas equações faltosas foram:
 e 
Essa escolha define uma spline natural. Fora do intervalo delimitado, a
spline é linear ou bastante próxima de uma função linear.
Relacionando todas as condições e fazendo as seguintes substituições:
As seguintes fórmulas para os coeficientes são obtidas:
 
 
Os valores de e de são constantes e nós conhecemos, mas ainda
é necessário descobrir todos os 
Pela continuidade das derivadas, monta-se o sistema linear AX=B tal
que:
A
Continuação do sistema linear:
Algoritmo
function ret = spline (x, y, xi)
 x = x(:);
 n = length (x);
 if (n < 3)
 error ("spline: Requer pelo menos 3 pontos!"); 
 endif
 ndy = ndims (y);
 szy = size (y);
 if (ndy == 2 && (szy(1) == 1 || szy(2) == 1))
 if (szy(1) == 1)
 a = y.';
 else
 a = y;
 szy = fliplr (szy);
 endif
 else
 a = reshape (y, [prod(szy(1:end-1)), szy(end)]).';
 endif
 complete = false;
 if (size (a, 1) == n + 2)
 complete = true;
 dfs = a(1,:);
 dfe = a(end,:);
 a = a(2:end-1,:);
 endif
 b = c = zeros (size (a));
 h = diff (x);
 idx = ones (columns (a), 1);
 if (complete)
 if (n == 3)
 dg = 1.5 * h(1) - 0.5 * h(2);
 c(2:n-1,:) = 1/dg(1);
 else
 dg = 2 * (h(1:n-2) .+ h(2:n-1));
 dg(1) = dg(1) - 0.5 * h(1);
 dg(n-2) = dg(n-2) - 0.5 * h(n-1);
 e = h(2:n-2);
 size(a)
 size(h)
 n
 g = 3 * diff (a(2:n,:)) ./ h(2:n-1,idx) ...
 - 3 * diff (a(1:n-1,:)) ./ h(1:n-2,idx);
 g(1,:) = 3 * (a(3,:) - a(2,:)) / h(2) ...
 - 3 / 2 * (3 * (a(2,:) - a(1,:)) / h(1) - dfs);
 g(n-2,:) = 3 / 2 * (3 * (a(n,:) - a(n-1,:)) / h(n-1) - dfe) ...
 - 3 * (a(n-1,:) - a(n-2,:)) / h(n-2);
 c(2:n-1,:) = spdiags ([[e(:); 0], dg, [0; e(:)]],
 [-1, 0, 1], n-2, n-2) \ g;
 endif
 c(1,:) = (3 / h(1) * (a(2,:) - a(1,:)) - 3 * dfs
 - c(2,:) * h(1)) / (2 * h(1));
 c(n,:) = - (3 / h(n-1) * (a(n,:) - a(n-1,:)) - 3 * dfe
+ c(n-1,:) * h(n-1)) / (2 * h(n-1));
 b(1:n-1,:) = diff (a) ./ h(1:n-1, idx) ...
- h(1:n-1,idx) / 3 .* (c(2:n,:) + 2 * c(1:n-1,:));
 d = diff (c) ./ (3 * h(1:n-1, idx));
 else
 g = zeros (n-2, columns (a));
 g(1,:) = 3 / (h(1) + h(2)) ...
* (a(3,:) - a(2,:) - h(2) / h(1) * (a(2,:) - a(1,:)));
 g(n-2,:) = 3 / (h(n-1) + h(n-2)) ...
* (h(n-2) / h(n-1) * (a(n,:) - a(n-1,:)) - (a(n-1,:) - a(n-2,:)));
 if (n > 4)
 g(2:n - 3,:) = 3 * diff (a(3:n-1,:)) ./ h(3:n-2,idx) ...
 - 3 * diff (a(2:n-2,:)) ./ h(2:n - 3,idx);
 dg = 2 * (h(1:n-2) .+ h(2:n-1));
 dg(1) = dg(1) - h(1);
 dg(n-2) = dg(n-2) - h(n-1);
 ldg = udg = h(2:n-2);
 udg(1) = udg(1) - h(1);
 ldg(n - 3) = ldg(n-3) - h(n-1);
 c(2:n-1,:) = spdiags ([[ldg(:); 0], dg, [0; udg(:)]],
 [-1, 0, 1], n-2, n-2) \ g;
 elseif (n == 4)
 dg = [h(1) + 2 * h(2); 2 * h(2) + h(3)];
 ldg = h(2) - h(3);
 udg = h(2) - h(1);
 c(2:n-1,:) = spdiags ([[ldg(:);0], dg, [0; udg(:)]],
 [-1, 0, 1], n-2, n-2) \ g;
 
 else # n == 3
 
 dg = h(1) + 2 * h(2);
 c(2:n-1,:) = g/dg(1);
 endif
 c(1,:) = c(2,:) + h(1) / h(2) * (c(2,:) - c(3,:));
 c(n,:) = c(n-1,:) + h(n-1) / h(n-2) * (c(n-1,:) - c(n-2,:));
 b = diff (a) ./ h(1:n-1, idx) ...
- h(1:n-1, idx) / 3 .* (c(2:n,:) + 2 * c(1:n-1,:));
 d = diff (c) ./ (3 * h(1:n-1, idx));
 endif
 d = d(1:n-1,:);
 c = c(1:n-1,:);
 b = b(1:n-1,:);
 a = a(1:n-1,:);
 coeffs = [d(:), c(:), b(:), a(:)];
 ret = mkpp (x, coeffs, szy(1:end-1));
 if (nargin == 3)
 ret = ppval (ret, xi);
 endif
endfunction