Ejercicio: 1Eva_2021PAOII_T1 Interpolación para perfil de terreno
a. Plantee y desarrolle un polinomio P3(d1) de grado 3, que describa el perfil del terreno en el intervalo [0,1300] de distancias a la primera antena d1.
Para plantear el ejercicio, se observan los punto en el intervalo [0,1300], asi como los tamaños de paso.
Como el polinomio es de grado 3, solo se requien 3 tramos o 4 muestras, el intervalo presenta 5. Para incluir el punto de inicio, fin y máximo, las muestras seleccionadas son las que se muestran en la tabla siguiente:
Δd1 | distancia d1 | Perfil de Terreno | DifDiv1 | DifDiv2 | DifDiv3 |
350 | 0 | 85 | DifDiv11 = 0.02857 | DifDiv12 = -6.122E-05 | DifDiv13 = -1.123E-05 |
350 | 350 | 95 | DifDiv21 = -0.01428 | DifDiv22 = -1.127E-05 | |
600 | 700 | 90 | DifDiv31 = -0.025 |
||
1300 | 75 |
Los tamaños de paso, Δd1, son diferentes, por lo que se pueden usar el método de diferencias divididas de Newton o el Método de Lagrange.
Usando diferencias divididas de Newton, se completa la tabla con las expresiones:
DifDiv_{11} = \frac{95-85}{350-0} = 0.02857143 DifDiv_{21} = \frac{90-95}{700-350} = -0.01428 DifDiv_{31} = \frac{75-90}{1300-700} =-0.025 DifDiv_{12} = \frac{-0.01428-0.02857}{700-0} = -6.122E-05 DifDiv_{22} = \frac{-0.025-(-0.01428)}{1300-350} = -1.127E-05 DifDiv_{13} = \frac{(-1.127E-05) -(6.122E-05)}{1300-0} = 1.123E-05con lo que se puede construir el polinomio
p_{3}(d) = 85 + 0.02857(d-0)+(-6.122E5)(d-0)(d-350)+ -(1.123E-5)(d-0)(d-350)(d-700)simplificando el polinomio se obtiene:
p_{3}(d) = 85.0+ 0.05941 d - 0.0001015 d^2 +(3.842E-8)d^3b. Calcule el error sobre el o los datos que no se usaron en el intervalo
El valor no usado que estaba en el intervalo es d=1000, que se evalúa en p3(d) y se compara con el valor muestreado para ver el error del modelo
p_{3}(1000) = 85.0+ 0.05941 (1000) - 0.0001015 (1000)^2 +(3.842E-8)(1000)^3 p_{3}(1000) = 81.26 Error = |81.26 - 80| = 1.26encontrando que el error es de 1.26m de altitud del terreno a los 1000m de distancia desde la referencia a la izquierda.
c. Desarrolle y justifique una propuesta para disminuir los errores encontrados en el literal anterior, sobre el mismo intervalo, es decir obtiene un nuevo polinomio (use algoritmo).
En caso de requerir una mayor precisión, se puede aumentar el grado del polinomio al incluir el punto d=1000 no usado en el paso anterior.
Usando el algoritmo con Python se obtiene:
Tabla Diferencia Dividida [['i ', 'xi ', 'fi ', 'F[1]', 'F[2]', 'F[3]', 'F[4]', 'F[5]']] [[ 0.0000e+00 0.0000e+00 8.5000e+01 2.8571e-02 -6.1224e-05 3.1920e-08 2.1666e-11 0.0000e+00] [ 1.0000e+00 3.5000e+02 9.5000e+01 -1.4286e-02 -2.9304e-05 6.0086e-08 0.0000e+00 0.0000e+00] [ 2.0000e+00 7.0000e+02 9.0000e+01 -3.3333e-02 2.7778e-05 0.0000e+00 0.0000e+00 0.0000e+00] [ 3.0000e+00 1.0000e+03 8.0000e+01 -1.6667e-02 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00] [ 4.0000e+00 1.3000e+03 7.5000e+01 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00]] dDividida: [ 2.8571e-02 -6.1224e-05 3.1920e-08 2.1666e-11 0.0000e+00] polinomio: 2.16658863275405e-11*x*(x - 1000.0)*(x - 700.0)*(x - 350.0) + 3.19204604918891e-8*x*(x - 700.0)*(x - 350.0) - 6.12244897959184e-5*x*(x - 350.0) + 0.0285714285714286*x + 85.0 polinomio simplificado: 2.16658863275405e-11*x**4 - 1.24946064795689e-8*x**3 - 6.6683650518237e-5*x**2 + 0.0525123706702654*x + 85.0
junto a la gráfica:
d. Escriba sus conclusiones y recomendaciones sobre los resultados obtenidos entre los dos polinomios.
Dado que el error se considera mínimo en el intervalo, y el coeficiente del polinomio para x4 es del orden de 10-11 se recomendaría usar el polinomio de grado 3.
Instrucciones en Python
# Polinomio interpolación # Diferencias Divididas de Newton # Tarea: Verificar tamaño de vectores, # verificar puntos equidistantes en x import numpy as np import sympy as sym import matplotlib.pyplot as plt # INGRESO , Datos de prueba xi = np.array([0.0,350,700,1000,1300]) fi = np.array([85.0,95,90,80,75]) # PROCEDIMIENTO # Tabla de Diferencias Divididas Avanzadas titulo = ['i ','xi ','fi '] n = len(xi) ki = np.arange(0,n,1) tabla = np.concatenate(([ki],[xi],[fi]),axis=0) tabla = np.transpose(tabla) # diferencias divididas vacia dfinita = np.zeros(shape=(n,n),dtype=float) tabla = np.concatenate((tabla,dfinita), axis=1) # Calcula tabla, inicia en columna 3 [n,m] = np.shape(tabla) diagonal = n-1 j = 3 while (j < m): # Añade título para cada columna titulo.append('F['+str(j-2)+']') # cada fila de columna i = 0 paso = j-2 # inicia en 1 while (i < diagonal): denominador = (xi[i+paso]-xi[i]) numerador = tabla[i+1,j-1]-tabla[i,j-1] tabla[i,j] = numerador/denominador i = i+1 diagonal = diagonal - 1 j = j+1 # POLINOMIO con diferencias Divididas # caso: puntos equidistantes en eje x dDividida = tabla[0,3:] n = len(dfinita) # expresión del polinomio con Sympy x = sym.Symbol('x') polinomio = fi[0] for j in range(1,n,1): factor = dDividida[j-1] termino = 1 for k in range(0,j,1): termino = termino*(x-xi[k]) polinomio = polinomio + termino*factor # simplifica multiplicando entre (x-xi) polisimple = polinomio.expand() # polinomio para evaluacion numérica px = sym.lambdify(x,polisimple) # Puntos para la gráfica muestras = 101 a = np.min(xi) b = np.max(xi) pxi = np.linspace(a,b,muestras) pfi = px(pxi) # SALIDA np.set_printoptions(precision = 4) print('Tabla Diferencia Dividida') print([titulo]) print(tabla) print('dDividida: ') print(dDividida) print('polinomio: ') print(polinomio) print('polinomio simplificado: ' ) print(polisimple) # Gráfica plt.plot(xi,fi,'o', label = 'Puntos') ##for i in range(0,n,1): ## plt.axvline(xi[i],ls='--', color='yellow') plt.plot(pxi,pfi, label = 'Polinomio') plt.legend() plt.xlabel('xi') plt.ylabel('fi') plt.title('Diferencias Divididas - Newton') plt.show()