Document 7519311

Download Report

Transcript Document 7519311

Diferenciación e Integración
numérica
Programación Numérica
Diferenciación
La diferenciación numérica puede calcularse usando la
definición de derivada
f '  x0   lim
h 0
f  x0  h   f x0 
h
Tomando una h pequeña. Si h > 0 se llama fórmula de
diferencia progresiva, si h < 0 se llama fórmula de diferencia
regresiva.
Error
x0= 2 ln(x0)= 0.693147181 f'(x0)= 0.5
h
0.1
0.01
0.001
0.0001
f(x0+h)
( f(x0+h) - f(x0) ) /h
0.741937345
0.487901642
0.698134722
0.498754151
0.693647056
0.499875042
0.693197179
0.4999875
|inc -f'(x0)|
-0.012098358
-0.001245849
-0.000124958
-1.24996E-05
Preguntas rápidas
Obtenga la derivada de las siguientes funciones en el punto
especificado utilizando Excel o Matlab. Compárelas con el
valor obtenido analíticamente.
1. f(x) = 3x sen(2x), x = p/6
2. f(x) = 5ln(x + 1) – x2/5, x = 1.2
Fórmulas de diferencias divididas
hacia adelante
Primera derivada
f  xi 1   f  xi 
f '  xi  
h
f '  xi  
 f  xi  2   4 f  xi 1   3 f  xi 
2h
Segunda derivada
f ' ' xi  
f  xi  2   2 f  xi 1   f xi 
h2
f ' '  xi  
 f  xi 3   4 f  xi  2   5 f  xi 1   2 f  xi 
h2
Tercera derivada
f xi 3   3 f xi  2   3 f xi 1   f xi 
f ' ' ' x  
i
f ' ' '  xi  
h3
 3 f xi  4   14 f xi 3   24 f xi  2   18 f  xi 1   5 f  xi 
2h 3
Fórmulas de diferencias divididas
centradas
Primera derivada
f  xi 1   f  xi 1 
f '  xi  
h
f '  xi  
 f  xi  2   8 f  xi 1   8 f  xi 1   f  xi  2 
12h
Segunda derivada
f  xi 1   2 f  xi   f  xi 1 
f ' '  xi  
h2
f ' ' xi  
 f xi  2   16 f xi 1   30 f xi   16xi 1   f  xi  2 
12h 2
Tercera derivada
f xi  2   2 f xi 1   2 f xi 1   f xi  2 
f ' ' ' x  
i
f ' ' '  xi  
2h 3
 f  xi 3   8 f xi  2   13 f  xi 1   13 f  xi 1   8 f  xi  2   f  xi 3 
8h 3
Fórmulas de diferencias divididas
hacia atrás
Primera derivada
f  xi   f  xi 1 
f '  xi  
h
f '  xi  
3 f  xi   4 f  xi 1   f  xi  2 
2h
Segunda derivada
f ' ' xi  
f  xi   2 f xi 1   f xi  2 
h2
f ' '  xi  
2 f  xi   5 f  xi 1   4 f  xi  2   f  xi 3 
h2
Tercera derivada
f xi   3 f  xi 1   3 f xi  2   f xi 3 
f ' ' ' x  
i
f ' ' '  xi  
h3
5 f xi   18 f xi 1   24 f  xi  2   14 f xi 3   3 f  xi  4 
2h 3
Ejemplo
f (x) = – 0.1x4 – 0.15x3 – 0.5x2 – 0.25x+1.2
x i-2
x i-1
xi
x i+1
x i+2
0
0.25
0.5
0.75
1
f'(x)
0.5
centrada
hacia adelante
hacia atrás
1.20000
1.10352
0.92500
0.63633
0.20000
-0.91250
O(h)
error
O(h^2)
error
-0.934375
2.4%
-0.9125
0.0%
-1.154688
26.5% -0.859375
5.8%
-0.714063
21.7% -0.878125
3.8%
Datos no espaciados
regularmente
Para derivar datos no espaciados regularmente se utiliza la
siguiente fórmula. Se requiere conocer la función en tres
puntos.
f ' x   f xi 1 
f xi 
2 x  xi  xi 1

xi 1  xi xi 1  xi 1 
2 x  xi 1  xi 1

xi  xi 1 xi  xi 1 
f xi 1 
2 x  xi 1  xi
xi 1  xi 1 xi 1  xi 
Ejemplo
El flujo de calor en la interfaz suelo-aire puede calcularse con la ley
de Faraday
dT
qz  0   kC
dz z 0
Donde q = flujo de calor, k = coeficiente de difusividad térmica
(3.5x10-7),  = la densidad del suelo (1800), C = calor específico
del suelo (840).
10
12 13.5
Aire
20  1.25  3.75
f ' 0  13.5

Suelo
0  1.250  3.75
1.25
20  0  3.75
20  0  1.25
12
 10
1.25  01.25  3.75 3.75  03.75  1.25
= –1.333
q = 70.56
3.75
Integración numérica
A los métodos de integración se les llama cuadratura numérica.
Seleccionaremos un conjunto de nodos [x0, ..., xn] del intervalo
[a, b].
Después integramos un polinomio interpolante de Lagrange
n
P x    f  xi Li  x 
i 0
Se obtiene:
Donde
b
n
a
i 0
 f x dx   a f x 
i
ai   Li x 
b
a
i
Regla del trapecio
Utilizando un polinomio interpolante lineal de Lagrange.
Px  
x  x1  f x   x  x0  f x 
x0  x1  0 x1  x0  1
b  x  x 
x  x0  f x  dx
1




f
x
dx

f
x

a
a  x0  x1  0 x1  x0  1 
x  x 
h
 1 0  f  x0   f x1    f x0   f  x1 
2
2
b
Donde h = x1 – x0 =
Esta fórmula vale cuando
f(x) tiene valores positivos.
P1
f
Da valores exactos para
polinomios de grado 1.
x0 = a
x1 = b
Pregunta rápida
Muestre que se cumple la regla del trapecio
 x  x1 


x  x0 
a f x dx  a  x0  x1  f x0   x1  x0  f x1  dx
x  x 
h
 1 0  f  x0   f x1    f x0   f  x1 
2
2
b
b
Regla se Simpson
La regla se Simpson se obtiene suponiendo el segundo polinomios
de Lagrange con los nodos x0 = a, x2 = b, x1 = a + h, h = (b – a)/2.
b   x  x  x  x 
x  x0 x  x2  f x   x  x0 x  x1  f x  dx
1
2




f
x
dx

f
x

a
a  x0  x1 x0  x2  0 x1  x0 x1  x2  1 x2  x0 x2  x1  2 
h
  f x0   4 f x1   f  x2 
3
b
Donde se han
despreciado los términos
de error.
La fórmula es exacta para
polinomios de hasta
tercer grado.
f
x0 = a
x1
P3
x2 = b
Comparación
Comparación entre el valor exacto, la regla del trapecio y
la regla de Simpson para diferentes funciones en el
intervalo [0 , 2].
f(x)
Valuación exacta
Trapecio
De Simpson
x^2
2.667
4.000
2.667
x^4
6.400
16.000
6.667
1/(x + 1)
1.099
1.333
1.111
sqrt(1 + x2)
2.958
3.236
2.964
sen x
1.416
0.909
1.425
exp(x)
6.389
8.389
6.421
Regla de Simpson 3/8
Ajustando polinomios de Lagrange de orden 3 usando cuatro
puntos se llega a la regla de Simpson de 3/8
I 
b
a
3h
f  x    f x0   3 f x1   3 f  x2   f  x3 
8
También puede expresarse por:

f  x0   3 f  x1   3 f  x2   f x3 
I   f  x   b  a 
a
8
b
Esta regla es útil cuando el número de puntos es impar.
Integración numérica compuesta
Integrando ex por Simpson en [0,4]

4
0
e x dx 


2 0
e  4e 2 ´e 4  56.76958
3
El error es: 53.59815 – 56.76958 = –3.17143
Separando en dos integrales:

4
0
2
4
e dx   e dx   e x dx
x
0

x
2
 
1 0
1
e  4e  e 2  e 2  4e 3  e 4
3
3
1
 e 0  4e  2e 2  4e 3  e 4
3
 53.86385




Dividiendo en 4 intervalos

4
0
1
2
3
4
e dx   e dx   e dx   e dx   e x dx
x



0
x
1
x
 
 
2
x
3

3
1
1 0
1
e  4e 2  e  e  4e 2  e 2
6
6
5
7
1
1
 e 2  4e 2  e3  e3  4e 2  e 4
6
6
3
5
7
1
1
 e 0  4e 2  2e  4e 2  2e 2  4e 2  2e3  4e 2  e 4
3
 53.61622


El error es: 53.59815 – 53.61622 = –0.01807

Regla compuesta de Simpson
Teorema. Sea f C4[a, b], n par, h = (b – a)/n, y xj = a + jh para
cada j = 0, 1, 2, ... n . La regla de Simpson para n subintervalos
puede escribirse como:

b
a
n / 2 1
n/2

h
f x dx   f a   2  f x2 j   4 f x2 j 1   f b
3
j 0
j 0

y= f(x)
x0 = a
x2
x2j-1
x2j
x2j+1
xn = b
Regla compuesta del trapecio
Teorema. Sea f C4[a, b], n par, h = (b – a)/n, y xj = a + jh para
cada j = 0, 1, 2, ... n . La regla del trapecio para n subintervalos
puede escribirse como:

b
a
n 1

h
f x dx   f a   2 f x j   f b 
2
j 1

y= f(x)
x0 = a x1
xj-1
xj
xn–1 xn = b
Regla compuesta del punto
medio
Teorema. Sea f C4[a, b], n par, h = (b – a)/(n+2), y xj = a +
(j+1)h para cada j = –1, 0, 1, 2, ... n+1. La regla de compuesta
del punto medio para n subintervalos puede escribirse como:
 f xdx  2h f x 
n/2
b
a
j 0
2j
y= f(x)
x0 = a x0
x1
xj-1
xj
xj+1
xn
xn+1 = b
Datos con espaciamiento
irregular
Si los datos están espaciados de forma irregular, como en el caso de datos
experimentales, la integración puede llevarse a cabo mediante la aplicación de la
regla del trapecio a cada subintervalo.
I  h1
f  x0   f  x1 
f  xn 1   f xn 
f  x1   f  x2 
 h2
 ...  hn
2
2
2
Donde hi = ancho del segmento i.
Ejemplo
Determinar la distancia recorrida para los datos
siguientes:
t min
1
2
3.25
4.5
6
7
8
9
9.5
10
V m/s
5
6
5.5
7
8.5
8
6
7
7
5
t = [1 2 3.25 4.5 6 7 8 9 9.5 10];
v = [5 6 5.5 7 8.5 8 6 7 7 5];
suma = 0;
for i=2:length(t)
suma = suma + (t(i)-t(i-1))*(v(i-1)+v(i))/2;
end
suma
ans = 60.3750
Algoritmos Regla del trapecio
Algoritmos para la regla del trapecio de uno solo segmento
function trap(h, f0, f1)
trap = h*(f0+f1)/2
end
Algoritmos para la regla del trapecio de múltiples segmentos
function trap(h, n, f)
sum = f0;
for i = 1, n–1
sum = sum + 2*fi
end
sum = sum + fn
trap = h*sum/2
end
Algoritmos Regla simple de
Simpson
Regla de Simpson de 1/3
function simp13(h, f0, f1, f2)
simp13 = 2*h*(f0+4*f1+f2)/6
end
Regla de Simpson de 3/8
function simp38(h, f0, f1, f2, f3)
simp38 = 3*h*(f0+3*f1+3*f2+f3)/8
end
Regla de Simpson 1/3 múltiple
Function simp13m(h, n, f)
sum = f0
for i = 1, n–2, 2
sum = sum+4*fi+2*fi+1
end
sum = sum+4fn-1+fn
simp13m = h*sum/3
end
Algoritmos Regla compuesta de Simpson
Regla de Simpson de número de segmentos pares o impares
function simpint(a, b, n, f)
h = (b-a)/n
if n=1 then
sum=trap()
else
m = n
odd = n/2-int(n/2)
if odd>0 and n>1 then
sum = sum + simp38(h,fn-3,fn-2,fn-1,fn)
m = n-3
end
if m>1 then
sum = sum + simp13m(h, m, f)
end
end
simpint = sum
end
Ejemplo Trapecio
Sea la siguiente función:
f (x) = 0.2 + 25x – 200x2 + 675x3 – 900x4 + 400x5
Integrada en el intervalo de a = 0 a b = 0.8 con trapecio:
Valor real I = 1.64053333
f (a) = 0.2000 f (b) = 0.2320
I = h (f (b) – f (a) )/2 0.17280000
error = 89.47%
Ejemplo Simpson 1/3
Sea la siguiente función:
f (x) = 0.2 + 25x – 200x2 + 675x3 – 900x4 + 400x5
Integrada en el intervalo de a = 0 a b = 0.8 con trapecio:
Valor real I = 1.64053333
f (a) = 0.2
f ((a+b)/2) = 2.456 f (b) = 0.232
I = 0.8 (0.2+4(2.456)+0.232)/6 = 1.36746667 error = 16.6%
Ejemplo Simpson 3/8
Sea la siguiente función:
f (x) = 0.2 + 25x – 200x2 + 675x3 – 900x4 + 400x5
Integrada en el intervalo de a = 0 a b = 0.8 con trapecio:
Valor real I = 1.64053333
f (0) = 0.2
f (0.26667) = 1.432724
f (0.5333) = 3.487177
f (0.8) = 2.232
I = 0.8 (0.2+3(1.432724+3.487177)+ 2.232 )/8 = 1.519170
error = 7.4%
Ejemplo Simpson 1/3 y Simpson 3/8
Sea la siguiente función:
f (x) = 0.2 + 25x – 200x2 + 675x3 – 900x4 + 400x5
Integrada en el intervalo de a = 0 a b = 0.8 con 5 segmentos,
con trapecio 2 primeros y Simpson los 3 últimos:
Valor real I = 1.64053333
f (0) = 0.2
f (0.16) = 1.29692
f (0.32) = 1.74339
f (0.48) = 3.18601
f (0.64) = 3.18193 f (0.8) = 0.23200
Simpson 1/3:
I1/3 = 0.32*(0.2 +4(1.29692)+ 1.74339 )/6 = 0.3803237
Simpson 3/8
I3/8 = 0.48 (1.74339 +3(3.18601 + 3.18193 )+ 2.232 )/8
= 1.264754
I = 1.645077
error = 0.28%
x
0.00
0.12
0.22
0.32
0.36
0.40
0.44
0.54
0.64
0.70
0.80
f(x)
0.200000
1.309729
1.305241
1.743393
2.074903
2.456000
2.842985
3.507297
3.181929
2.363000
0.232000