Ejercicio: 1Eva2015TI_T2 Salida cardiaca
Solución presentada como introducción al tema de interpolación y solución de sistemas de ecuaciones. Como ejemplo, se usa el método de Lagrange.
# Gráfica de datos experimentales:
ti = [2,6,9,12,15,18]
yi = [0,1.5,3.2,4.1,3.4,2.0]literal a
Para seleccionar los puntos del polinomio de grado 2, se requiere observar la gráfica.
Se requiere desarrollar la parte analítica en papel y lápiz para al menos una de las gráficas. (tarea)
Donde se observa que se podría aproximar a una parábola usando los puntos extremos y alrededor de t=12.
ti = [2 ,12, 18]
yi = [0,4.1,2.0]con lo que se obtiene la gráfica siguiente:
Polinomio de Lagrange:
-0.0475*x**2 + 1.075*x - 1.96
2
- 0.0475⋅x + 1.075⋅x - 1.96
En la primera parte de la gráfica en el intervalo [2,12] presenta un mayor error respecto a los puntos experimentales.
Se intenta aumentando el punto alrededor de 6
ti = [2,6,12,18]
yi = [0,1.5,4.1,2.0]con un mejor resultado observado a simple vista
Polinomio de Lagrange:
-0.004444*x**3 + 0.094722*x**2 - 0.151667*x - 0.04
3 2
- 0.004444⋅x + 0.094722⋅x - 0.151667⋅x - 0.04
literal b
Usando el algoritmo se puede incorporar todos los puntos
Algoritmo en Python
# 1Eva2015TI_T2 Salida cardiaca
# Interpolacion de Lagrange
# divisoresL solo para mostrar valores denominador
import numpy as np
import sympy as sym
def interpola_Lagrange(xi,fi,vertabla=False,
precision=6, casicero = 1e-15):
'''
Interpolación con método de Lagrange
resultado: polinomio en forma simbólica
'''
# Matrices como arreglo, numeros reales
xi = np.array(xi,dtype=float)
fi = np.array(fi,dtype=float)
if vertabla==True:
print('Interpola con Lagrange')
# Polinomio con Lagrange
n = len(xi)
x = sym.Symbol('x')
polinomio = sym.S.Zero
divisorL = np.zeros(n, dtype = float)
for i in range(0,n,1):
termino = 1 # Termino de Lagrange
numerador = 1
denominador = 1
for j in range(0,n,1):
if (j!=i):
numerador = numerador*(x-xi[j])
denominador = denominador*(sym.UnevaluatedExpr(xi[i])-xi[j])
terminoLi = numerador/denominador
if vertabla==True:
txt0='' ; txt1='('; txt2=')'
if i>0:
txt0='+'
if n>2:
txt1=''; txt2=''
print(txt0,fi[i],'*'+txt1,numerador,txt2+'/'+ txt1,
denominador,txt2)
polinomio = polinomio + (terminoLi.doit())*fi[i]
denominador = denominador.doit()
divisorL[i] = denominador
# simplifica el polinomio
polisimple = polinomio.expand()
try: # intenta redondear coeficientes a precisiion
polisimple = redondea_coef(polisimple, precision,casicero)
except NameError:
print("redondea_coef() no se encuentra definida.")
if vertabla==True:
print('polinomio simplificado')
print(polisimple)
return(polisimple)
def redondea_coef(ecuacion, precision=6,casicero = 1e-15):
''' redondea coeficientes de términos suma de una ecuacion
'''
tipo = type(ecuacion)
tipo_eq = False
if tipo == sym.core.relational.Equality:
RHS = ecuacion.rhs
ecuacion = ecuacion.lhs
tipo = type(ecuacion)
tipo_eq = True
if tipo == sym.core.add.Add: # términos suma de ecuacion
term_sum = sym.Add.make_args(ecuacion)
ecuacion = sym.S.Zero
for term_k in term_sum:
# factor multiplicativo de termino suma
term_mul = sym.Mul.make_args(term_k)
producto = sym.S.One
for factor in term_mul:
if not(factor.has(sym.Symbol)):
factor = np.around(float(factor),precision)
if (abs(factor)%1)<casicero: # si es entero
factor = int(factor)
producto = producto*factor
ecuacion = ecuacion + producto
if tipo == sym.core.mul.Mul: # termino único, busca factores
term_mul = sym.Mul.make_args(ecuacion)
producto = sym.S.One
for factor in term_mul:
if not(factor.has(sym.Symbol)):
factor = np.around(float(factor),precision)
if (abs(factor)%1)<casicero: # si es entero
factor = int(factor)
producto = producto*factor
ecuacion = producto
if tipo == float: # si es entero
if (abs(ecuacion)%1)<casicero:
ecuacion = int(ecuacion)
if tipo_eq:
ecuacion = sym.Eq(ecuacion,RHS)
return(ecuacion)
# PROGRAMA DE PRUEBA -----------------------------
# INGRESO
# Gráfica de datos experimentales:
ti = [2,6,9,12,15,18]
yi = [0,1.5,3.2,4.1,3.4,2.0]
xi = [2,6,12,18]
fi = [0,1.5,4.1,2.0]
# PROCEDIMIENTO
polisimple = interpola_Lagrange(xi,fi,vertabla=True)
x = sym.Symbol('x')
px = sym.lambdify(x,polisimple)
# SALIDA
print('Polinomio de Lagrange: ')
print(polisimple)
print()
sym.pprint(polisimple)
# Gráfica --------------
import matplotlib.pyplot as plt
muestras = 21 # resolución gráfica
a = np.min(ti) # intervalo [a,b]
b = np.max(ti)
xk = np.linspace(a,b,muestras)
yk = px(xk)
plt.plot(ti,yi,'o', label = '[xi,fi]')
plt.plot(xk,yk, label = 'p(x)')
plt.legend()
plt.xlabel('xi')
plt.ylabel('fi')
plt.title('Interpolación Lagrange')
plt.show()
Observe que si solamente usa los 3 primeros puntos, la curva obtenida aumenta el error para una parte del intervalo.
literal c
No se desarrolla el tema de integrales. (tarea)
Revisar el tema detallado de integrales, considerando todos los puntos experimentales.