s1Eva_IIT2007_T2 Aplicar Gauss-Seidel

Probando solución con Jacobi, enviado como tarea

Desarrrollo Analítico

– Verificar que la matriz es diagonal dominante

No es necesario realizar el pivoteo por filas, ya la matriz tiene la diagonal dominante.

A = [[7.63, 0.30, 0.15,  0.50, 0.34, 0.84],
     [0.38, 6.40, 0.70,  0.90, 0.29, 0.57],
     [0.83, 0.19, 8.33,  0.82, 0.34, 0.37],
     [0.50, 0.68, 0.86, 10.21, 0.53, 0.70],
     [0.71, 0.30, 0.85,  0.82, 5.95, 0.55],
     [0.43, 0.54, 0.59,  0.66, 0.31, 9.25]]

B = [ -9.44, 25.27, -48.01, 19.76, -23.63, 62.59]

– revisar el número de condición

cond(A) = || A||.||A-1||

El número de condición no es «muy alto»,  los valores de la diagonal son los mayores en toda la fila, por lo que el sistema converge.

>>> np.linalg.cond(A)
2.0451853966291011

– realizar iteraciones

Dado que no se establece en el enunciado el vector inicial, se usará el vector cero. La tolerancia requerida es 10-5
X0 = [ 0. 0. 0. 0. 0. 0.]

iteración 1

x0 = (-9.44 -0.30(0) -0.15(0) 
      -0.50(0) -0.34(0) -0.84(0))/7.63 
   = -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(0) -0.70(0) 
      -0.90(0) -0.29(0) -0.57(0))/6.40 
   = 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(0) -0.19(0) -0.82(0)
      -0.34(0) -0.37(0))/8.33 
   = -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(0) -0.68(0) -0.86(0)
      -0.53(0) -0.70(0))/10.21 
   = 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(0) -0.30(0) -0.85(0)
      -0.82(0) -0.55(0))/5.95 
   = -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(0) -0.54(0) -0.59(0)
      -0.66(0) -0.31(0))/9.25 
   = 62.59/9.25 = 6.76648649
X1 = [-1.23722149  3.9484375  -5.7635054   1.93535749 -3.97142857  6.76648649]
diferencia = X1 - X0 = X1
errado = max(|diferencia|) = 6.76648649

iteración 2

x0 = (-9.44 -0.30(3.9484375) -0.15(-5.7635054) 
      -0.50(1.93535749) -0.34(-3.97142857) -0.84(6.76648649))/7.63 
   = -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(-1.23722149) -0.70(-5.7635054) 
      -0.90(1.93535749) -0.29(-3.97142857) -0.57(6.76648649))/6.40 
   = 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(-1.23722149) -0.19(3.9484375)
      -0.82(1.93535749) -0.34(-3.97142857) -0.37(6.76648649))/8.33 
   = -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(-1.23722149) -0.68(3.9484375)
      -0.86(-5.7635054) -0.53(-3.97142857) -0.70(6.76648649))/10.21 
   = 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(-1.23722149) -0.30(3.9484375) 
      -0.85(-5.7635054) -0.82(1.93535749) -0.55(6.76648649))/5.95 
   = -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(-1.23722149) -0.54(3.9484375) 
      -0.59(-5.7635054) -0.66(1.93535749) -0.31(-3.97142857))/9.25 
   = 62.59/9.25 = 6.76648649
X1 = [-1.97395113  3.95743644 -6.05925771  1.96068604 -4.09171178  6.95612152]
diferencia = X1 - X0 = [-0.73672964,  0.00899894, -0.29575231,  0.02532855, -0.12028321, 0.18963504]
errado = max(|diferencia|) = 0.736729635697

iteración 3

x0 = (-9.44 -0.30(3.95743644) -0.15(-6.05925771) 
      -0.50(1.96068604) -0.34(-4.09171178) -0.84(6.95612152))/7.63 
   = -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(-1.97395113) -0.70(-6.05925771) 
      -0.90(1.96068604) -0.29(-4.09171178) -0.57(6.95612152))/6.40 
   = 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(-1.97395113) -0.19(3.95743644) 
      -0.82(1.96068604) -0.34(-4.09171178) -0.37(6.95612152))/8.33 
   = -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(-1.97395113) -0.68(3.95743644) 
     -0.86(-6.05925771) -0.53(-4.09171178) -0.70(6.95612152))/10.21 
   = 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(-1.97395113) -0.30(3.95743644) 
      -0.85(-6.05925771) -0.82(1.96068604) -0.55(6.95612152))/5.95 
   = -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(-1.97395113) -0.54(3.95743644) 
      -0.59(-6.059257710) -0.66(1.96068604) -0.31(-4.09171178))/9.25 
   = 62.59/9.25 = 6.76648649
X1 = [-1.98566781  4.0185268  -5.9920623   2.01431955 -3.98302285  7.01093224]
diferencia = X1 - X0 = [-0.01171668,  0.06109037,  0.0671954 ,  0.05363351,  0.10868893, 0.05481072]
errado = max(|diferencia|) = 0.108688931048

Desarrollo numérico con Python

Se verifica el resultado obtenido realizando A.Xi y comparando con el vecto B
en la tabla se usa el signo de errado para la gráfica.

X0:  [ 0.  0.  0.  0.  0.  0.]
Xi:  [-1.23722149  3.9484375  -5.7635054   1.93535749 -3.97142857  6.76648649]
errado:  6.76648648649
Xi:  [-1.97395113  3.95743644 -6.05925771  1.96068604 -4.09171178  6.95612152]
errado:  -0.736729635697
Xi:  [-1.98566781  4.0185268  -5.9920623   2.01431955 -3.98302285  7.01093224]
errado:  0.108688931048
Xi:  [-2.00378293  3.99452422 -6.00443878  1.99576482 -4.0067623   6.9961552 ]
errado:  -0.0240025781157
Xi:  [-1.99869529  4.00195452 -5.99863447  2.00153846 -3.99769932  7.00130746]
errado:  0.00906298532211
Xi:  [-2.00045097  3.99933614 -6.00047801  1.99948184 -4.00078219  6.99955127]
errado:  -0.00308287453133
Xi:  [-1.99984629  4.00022733 -5.99983706  2.00017793 -3.99973154  7.00015339]
errado:  0.0010506528253
Xi:  [-2.00005265  3.9999222  -6.00005579  1.99993915 -4.00009178  6.9999475 ]
errado:  -0.0003602430897
Xi:  [-1.99998199  4.00002662 -5.99998091  2.00002082 -3.99996859  7.00001796]
errado:  0.000123195668206
Xi:  [-2.00000616  3.99999089 -6.00000653  1.99999287 -4.00001075  6.99999385]
errado:  -4.21623029112e-05

respuesta de A.X=B : 
[-2.00000616  3.99999089 -6.00000653  1.99999287 -4.00001075  6.99999385]
iteraciones:  10

A.Xi:  [ -9.44006312  25.26992175 -48.01007303  19.75990236 -23.63008584
  62.58992368]
   B: [ -9.44 25.27 -48.01 19.76 -23.63 62.59]

el gráfico de los errores vs iteraciones es:

s1Eva_IT2011_T1 Encontrar α en integral

Desarrollo Analítico

Se iguala la ecuación al valor buscado = 10, y se resuelve

\int_{\alpha}^{2\alpha} x e^{x}dx = 10

siendo: μ = x , δv = ex, δu = δx , v = ex

\int u dv = uv - \int v \delta u xe^x \Big|_{\alpha}^{2 \alpha} - \int_{\alpha}^{2\alpha} e^{x}dx - 10 = 0 2\alpha e^{2 \alpha} -\alpha e^{\alpha} - (e^{2\alpha} - e^{\alpha}) - 10 = 0 (2\alpha-1)e^{2 \alpha}+ (1-\alpha) e^{\alpha} - 10 = 0

la función a usar en el método es

f(\alpha) = (2\alpha-1)e^{2 \alpha}+ (1-\alpha)e^{\alpha} -10

Se obtiene la derivada para el método de Newton Raphson

f'(\alpha) = 2e^{2 \alpha} + 2(2\alpha-1)e^{2 \alpha} - e^{\alpha} + (1-\alpha) e^{\alpha} f'(\alpha) = (2 + 2(2\alpha-1))e^{2 \alpha} +(-1 + (1-\alpha)) e^{\alpha} f'(\alpha) = 4\alpha e^{2 \alpha} -\alpha e^{\alpha}

la fórmula para el método de Newton-Raphson

\alpha_{i+1} = \alpha_i - \frac{f(\alpha)}{f'(\alpha)}

se prueba con α0 = 1, se evita el valor de cero por la indeterminación que se da por f'(0) = 0

iteración 1:

f(1) = (2(1)-1)e^{2(1)}+ (1-(1))e^{(1)} -10 f'(1) = 4(1) e^{2 (1)} -(1) e^{(1)} \alpha_{1} = 1- \frac{f(1)}{f'(1)}

\alpha_{1} = 1.0973
error = 0.0973

iteración 2:
f(2) = (2(1.0973)-1)e^{2(1.0973)}+ (1-(1.0973))e^{(1.0973)} -10

f'(2) = 4(1.0973) e^{2 (1.0973)} -(1.0973) e^{(1.0973)} \alpha_{2} = 1.0973 - \frac{f(1.0973)}{f'(1.0973)}

\alpha_{2} = 1.0853
error = 0.011941

iteración 3:
f(3) = (2(1.0853)-1)e^{2(1.0853)}+ (1-(1.0853))e^{(1.0853)} -10

f'(3) = 4(1.0853) e^{2 (1.0853)} -(1.0853) e^{(1.0853)} \alpha_{3} = 1.0853- \frac{f(1.0853)}{f'(1.0853)}

\alpha_{3} = 1.0851
error = 0.00021951

[  xi,          xnuevo,      f(xi),      f'(xi),       tramo ]
[  1.           1.0973      -2.6109      26.8379       0.0973]
[  1.0973e+00   1.0853e+00   4.3118e-01   3.6110e+01   1.1941e-02]
[  1.0853e+00   1.0851e+00   7.6468e-03   3.4836e+01   2.1951e-04]
[  1.0851e+00   1.0851e+00   2.5287e-06   3.4813e+01   7.2637e-08]
raiz:  1.08512526549

se obtiene el valor de la raíz con 4 iteraciones, con error de aproximación de 7.2637e-08

Desarrollo con Python

s1Eva_IT2009_T1 Demanda de producto

Desarrollo analítico

– igualar la ecuación al valor buscado, 80

200 t e^{-0.75t} = 80

– forma estándar de la ecuación para el método f(x) = 0:

f(t) = 200 t e^{-0.75t} - 80

– derivada de la ecuación

f'(t) = 200 e^{-0.75t} + 200 t (-0.75) e^{-0.75t} f'(t) = 200 e^{-0.75t}(1-0.75t)

– fórmula del método de Newton-Raphson

t_{i+1} = t_i - \frac{f(t)}{f'(t)}

– Punto inicial. Como la variable es t, tiempo, el rango de análisis es t>0
El valor inicial de búsqueda se selecciona t0 = 1

iteración 1:

t_{1} = 1 - \frac{200 (1) e^{-0.75(1)} - 80}{200 e^{-0.75(1)}(1-0.75(1))}

error = 0.6128

iteración 2:

t_{2} = 0.3872 - \frac{200 (0.3872) e^{-0.75(0.3872)} - 80}{200 e^{-0.75(0.3872)}(1-0.75(0.3872))}

error = 0.208

iteracion 3:

t_{3} = 0.5952 - \frac{200 (0.5952) e^{-0.75(0.5952)} - 80}{200 e^{-0.75(0.5952)}(1-0.75(0.5952))}

error = 5.3972e-02

tabla de iteraciones

[  xi,       xnuevo,   f(xi),   f'(xi),    tramo]
[  1.        0.3872   14.4733   23.6183   0.6128]
[  0.3872    0.5952  -22.0776  106.1511   0.208 ]
[  5.9518e-01   6.4916e-01  -3.8242e+00   7.0855e+01   5.3972e-02]
[  6.4916e-01   6.5252e-01  -2.1246e-01   6.3069e+01   3.3686e-03]
[  6.5252e-01   6.5254e-01  -7.9031e-04   6.2600e+01   1.2625e-05]
raiz:  0.65253630273

Se obtiene el valor de la raíz con 5 iteraciones, y un error de 1.2625e-05

Desarrollo con Python:

Tarea: desarrollar el literal b

s1Eva_IIT2011_T2 Sistema de Ecuaciones, diagonal dominante

1. Desarrollo analítico

1.1 Solución iterativa usando el médodo de Jacobi

El sistema de ecuaciones

\begin{cases} -2x+5y+9z=1\\7x+y+z=6\\-3x+7y-z=-26\end{cases}

cambia a su forma matricial AX=B

\begin{bmatrix} -2 & 5 & 9 \\ 7 & 1 & 1 \\ -3 & 7 &-1 \end{bmatrix} \begin{bmatrix} x_0 \\ x_1 \\ x_2 \end{bmatrix} = \begin{bmatrix} 1 \\ 6 \\ -26 \end{bmatrix}

para usarla en el algoritmo se intercambian filas, pivoteo, buscando hacerla diagonal dominante:

\begin{bmatrix} 7 & 1 & 1 \\ -3 & 7 &-1\\ -2 & 5 & 9 \end{bmatrix} \begin{bmatrix} x_0 \\ x_1 \\ x_2 \end{bmatrix} = \begin{bmatrix} 6 \\ -26 \\1 \end{bmatrix}

las ecuaciones para el algoritmo se obtienen despejando una variable diferente en cada ecuación.

\begin{cases}x=(6-y-z)/(7)\\y=(-26+3x+z)/(7)\\z=(1+2x-5y)/(9)\end{cases}

Dado el vector inicial X(0) = [0, 0, 0], y una tolerancia de 0.0001, e desarrollan al menos 3 iteraciones:

Nota. El super índice entre paréntesis denota el número de la iteración. El subíndice denota el elemento del vector. El error se verifica para éste ejercicio como el mayor valor de la diferencia entre iteraciones consecutivas. Analogía al video del acople entre las aeronaves. Si una coordenada aún no es la correcta …. menor a lo tolerado, pues NO hay acople…

Iteración 1

\begin{cases}x=(6-(0)-(0))/(7)\\y=(-26+3(0)+(0))/(7)\\z=(1+2(0)-5(0))/(9)\end{cases}
X(1) = [6/7, -26/7, 1/9] 
    = [0.8571, -3,7142, 0.1111]
diferencia(1) = |X(1) - X(0)|
             = |[0.8571,-3,7142,0.1111] - [0,0,0]|
             = |[0.8571,-3,7142,0.1111]|
errado(1) = maximo|[0.8571,-3,7142,0.1111]|
         = 26/7 ≅ 3.71, 
error es más alto que lo tolerado

Iteración 2

\begin{cases}x=(6-(-26/7)-(1/9))/(7)\\y=(-26+3(6/7)+(1/9))/(7)\\z=(1+2(6/7)-5(-26/7))/(9)\end{cases}
X(2) = [1.3718, -3.3310 ,  2,3650]
diferencia(2) = |X(2) - X(1)|
             = |[1.3718, -3.3310 , 2,3650] - [0.8571, -3,7142, 0.1111]| 
             = |[0.5146, 0.3832, 2.2539]|
errado(2) ≅ 2.2538, 
el error disminuye, pero es más alto que lo tolerado

Iteración 3

\begin{cases}x=(6-(-3.3310)-( 2,3650))/(7)\\y=(-26+3x+z)/(7)\\z=(1+2x-5y)/(9)\end{cases}
X(3) = [0.9951, -2.7884, 2.2665]
diferencia(3) = |[0.9951, -2.7884, 2.2665] - [1.3718, -3.3310 ,  2,3650]|
             = |[-0.3767, 0.5425, -0.0985]|
errado(3) ≅ 0.5425, 
el error disminuye, pero es más alto que lo tolerado

Observación: Si el error disminuye en cada iteración, se puede intuir que se converge a la solución. Se puede continuar con la 4ta iteración…


2. Solución numérica usando Python

Usando el algoritmo desarrollado en clase se obtiene la respuesta con 13 iteraciones:

Xi, errado
[ 0.85714286 -3.71428571  0.11111111] 3.71428571429
[ 1.37188209 -3.33106576  2.36507937] 2.25396825397
[ 0.99514091 -2.78846777  2.26656589] 0.542597991578
[ 0.93170027 -2.96400162  1.8814023 ] 0.385163589245
[ 1.0117999  -3.04621384  1.96482318] 0.0834208883016
[ 1.01162724 -2.99996816  2.02829656] 0.0634733731283
[ 0.99595309 -2.99097453  2.00256614] 0.0257304176216
[ 0.99834406 -3.0013678   1.99408654] 0.0103932672897
[ 1.00104018 -3.00155447  2.0003919 ] 0.00630536414945
[ 1.00016608 -2.99949822  2.00109475] 0.00205624814582
[ 0.99977193 -2.99977243  1.99975814] 0.00133660433338
[ 1.00000204 -3.0001323   1.99982289] 0.000359867509479
[ 1.0000442  -3.00002443  2.00007395] 0.000251063280064
[ 0.99999292 -2.99997049  2.00002339] 5.39347670632e-05
iteraciones:  13

La gráfica muestra las coordenadas de las aproximaciones, observe el espiral que se va cerrando desde el punto inicial en verde al punto final en rojo.

Tarea: convierta el algoritmo a una función, así es más sencillo usar en otros ejercicios.
Se debe controlar el número de iteraciones para verificar la convergencia con iteramax.
NO es necesario almacenar los valores de los puntos, solo el último valor y el contador itera permite determinar si el sistema converge.

# 1ra Evaluación II Término 2011
# Tema 2. Sistema ecuaciones Jacobi
import numpy as np

# INGRESO
A = np.array([[ 7.0, 1, 1],
              [-3.0, 7,-1],
              [-2.0, 5, 9]])

B = np.array([6.0, -26.0, 1.0])

X = np.array([0.0, 0.0,  0.0],dtype=float)
tolera = 0.0001
iteramax = 100

# PROCEDIMIENTO
tamano = np.shape(A)
n = tamano[0]
m = tamano[1]
Xi1 = np.zeros(n, dtype=float)
errado = 2*tolera
itera = 0
while not(errado<=tolera or itera>=iteramax):
    i = 0
    while not(i>=n):
        j = 0
        nuevo = B [i]
        while not(j>=m):
            if (i!=j):
                nuevo = nuevo - A[i,j]*X[j]
            j = j+1
        Xi1[i] = nuevo/A[i,i]
        i = i+1
    diferencia = np.abs(X - Xi1)
    errado = np.max(diferencia)
    X = np.copy(Xi1)
    itera = itera +1
    print(Xi1, errado)

# SALIDA
print('iteraciones: ', itera-1)

Tarea: Realice las modificaciones necesarias para usar el algoritmo de Gauss-Seidel. Luego compare resultados.


Revisión de resultados

Si se usan las ecuaciones sin cambiar a diagonal dominante, como fueron presentadas en el enunciado, el algoritmo Jacobi NO converge, el error aumenta, y muy rápido como para observar la espiral hacia afuera en una gráfica.

Xi, errado
[ -0.5   6.   26. ] 26.0
[ 131.5  -16.5   69.5] 132.0
[ 271. -984. -484.] 967.5
[-4638.5 -1407.  -7675. ] 7191.0
[-38055.5  40150.5   4092.5] 41557.5
[ 118792.  262302.  395246.] 391153.5
[ 2434361.5 -1226784.   1479764. ] 2315569.5
[  3591977.5 -18520288.5 -15890546.5] 17370310.5
....

s1Eva_IT2018_T1 Tanque esférico canchas deportivas

a) Para seleccionar el rango para h=[a,b], se observa que el tanque puede estar vacio, a=0 o lleno al máximo, b=2R = 2(3)=6, con lo que se obtiene:

h =[0.0, 6.0]

conociendo la proporción con el valor máximo, se tiene un valor inicial para h0 para las iteraciones.

Vmax = \frac{\pi}{3} (2R)^2 (3R-2R) Vmax = \frac{4\pi }{3}R^{3}= 113.09 h_0 = (6)*30/113.09 = 1.59

b) Usar Newton-Raphson con tolerancia 1e-6

f(h) = \frac{\pi }{3}h^2 (3(3)-h)-30 f(h) = \frac{\pi }{3}(9h^2 -h^3-28.647) f'(h) = \frac{\pi }{3}(18h-3h^2)

el punto siguiente de iteración es:

h_{i+1} = h_{i} -\frac{f(h)}{f'(h)} = h_{i}-\frac{ \frac{\pi }{3}(9h^2 -h^3-28.647)}{ \frac{\pi }{3}(18h-3h^2)} h_{i+1} = h_{i} -\frac{(9h^2 -h^3-28.647)}{(18h-3h^2)}

con lo que se realiza la tabla de iteraciones:

hi hi+1 error orden
1.590 2.061 0.47 10-1
2.061 2.027 -0.034 10-2
2.027 2.02686 -0.00014 10-4
2.02686 2.0268689 -2.32E-09 10-9

En valor práctico 2.028 m usando flexómetro, a menos que use medidor laser con precisión 10-6 usará más dígitos con un tanque de 6 metros de altura ganará una precisión de una gota de agua para usar en duchas o regar el césped .

c) El orden de convergencia del error observando las tres primeras iteraciones es cuadrático

Tarea: Realizar los cálculos con Python, luego aplique otro método. Añada sus observaciones y conclusiones.

s1Eva_IT2018_T3 Temperatura en nodos de placa

a) Plantear el sistema de ecuaciones. Usando el promedio para cada nodo interior:

a=\frac{50+c+100+b}{4} b=\frac{a+30+50+d}{4} c=\frac{a+60+100+d}{4} d=\frac{b+60+c+30}{4}

que reordenando se convierte en:

4a=150+c+b 4b=a+80+d 4c=a+160+d 4d=b+c+90

simplificando:

4a-b-c= 150 a-4b+d = -80 a-4c+d = -160 b+c-4d = -90

que a forma matricial se convierte en:

A = [[ 4, -1, -1, 0.0],
     [ 1, -4,  0, 1.0],
     [ 1,  0, -4, 1.0],
     [ 0,  1,  1,-4.0]]
B = [[ 150.0],
     [ -80.0],
     [-160.0],
     [ -90.0]]

Observación: la matriz A ya es diagonal dominante, no requiere pivotear por filas.  Se aumentó el punto decimal a los valores de la matriz A y el vectot B  para que sean considerados como números reales.

El número de condición es: np.linalg.cond(A) = 3.0

que es cercano a 1 en un orden de magnitud, por lo que la solución matricial es «estable» y los cambios en los coeficientes afectan proporcionalmente a los resultados. Se puede aplicar métodos iterativos sin mayores inconvenientes.

b y c) método de Jacobi para sistemas de ecuaciones, con vector inicial

 
X(0) = [[60.0],
        [40], 
        [70],
        [50]] 

reemplazando los valores iniciales en cada ecuación sin cambios.

iteración 1
a=\frac{50+70+100+40}{4} = 65

b=\frac{60+30+50+50}{4} = 47.5 c=\frac{60+60+100+50}{4} = 67.5 d=\frac{40+60+70+30}{4} = 50
X(1) = [[65],
        [47.5],
        [67.5],
        [50]]

vector de error = 
     [|65-60|,
      |47.5-40|,
      |67.5-70|,
      |50-50|]
  = [|5|,
     |7.5|,
     |-2.5|,
     |0|]
errormax = 7.5

iteración 2
a=\frac{50+67.5+100+47.5}{4} = 66.25

b=\frac{65+30+50+50}{4} = 48.75 c=\frac{65+60+100+50}{4} = 68.75 d=\frac{47.5+60+67.7+30}{4} = 51.3
X(2) = [[66.25],
        [48.75],
        [68.75],
        [51.3]]

vector de error = 
       [|66.25-65|,
        |48.75-47.5|,
        |68.75-67.5|,
           |51.3-50|] 
  = [|1.25|,
     |1.25|,
     |1.25|,
     |1.3|]
errormax = 1.3

iteración 3
a=\frac{50+68.75+100+48.75}{4} = 66.875

b=\frac{66.25+30+50+51.3}{4} = 49.38 c=\frac{66.25+60+100+51.3}{4} = 69.3875 d=\frac{48.75+60+68.75+30}{4} = 51.875
X(2) = [[66.875],
        [49.38],
        [69.3875],
        [51.875]]

vector de error = 
      [|66.875-66.25|,
       |49.38-48.75|,
       |69.3875-68.75|,
       |51.875-51.3|]
    = [|0.655|,
       |0,63|,
       |0.6375|,
        |0.575|]
errormax = 0.655
con error relativo de:
100*0.655/66.875 = 0.97%

siguiendo las iteraciones se debería llegar a:

>>> np.linalg.solve(A,B)
array([[ 67.5],
       [ 50. ],
       [ 70. ],
       [ 52.5]])

s1Eva_IT2018_T4 El gol imposible

Tabla de datos:

ti = [0.00, 0.15, 0.30, 0.45, 0.60, 0.75, 0.90, 1.05, 1.20]
xi = [0.00, 0.50, 1.00, 1.50, 1.80, 2.00, 1.90, 1.10, 0.30]
yi = [0.00, 4.44, 8.88,13.31,17.75,22.19,26.63,31.06,35.50]
zi = [0.00, 0.81, 1.40, 1.77, 1.91, 1.84, 1.55, 1.03, 0.30]

Observe que, un gol simplificado con física básica es un tiro parabólico cuya trayectoria se compone de movimientos en los ejes, Y y Z. Sin embargo, lo «imposible» del gol mostrado es añadir el movimiento en X. (Referencia de la imagen en el enunciado)

El movimiento de «profundidad» o dirección hacia el arco y(t) es semejante a un polinomio de primer grado, y el movimiento de «altura» z(t) es un polinomio de segundo grado. El movimiento «imposible» en el eje X, podría ser descrito por un polinomio de segundo o mayor grado.

a) Encontrar t para altura máxima, que se encuentra al igualar la derivada dz/dt a cero. Para interpolar el polinomio z(t), de segundo grado, se puede usar tres puntos de los sugeridos: 0, 0.3 y 0.6, que además son equidistantes en t (facilita usar cualquier método de interpolación).

Por ejemplo, con diferencias finitas avanzadas:

t z(t) d1 d2 d3
0.00 0.00 1.40 -0.89
0.30 1.40 0.51
0.60 1.91
z(t) = 0 + 1.40\frac{(t-0)}{1!(0.3)} - 0.89 \frac{(t-0)(t-0.3)}{2!(0.3)^2} = 0 + 4.66 t - 4.94(t^2-0.3t) = 4.66 t + 1.48 t - 4.94 t^2 z(t) = 6.42 t - 4.94 t^2

para encontrar el valor máximo se encuentra \frac{dz}{dt} = 0

\frac{dz}{dt} = 6.42 - 2(4.94) t 6.42 - 2(4.94) t = 0 t = \frac{6.42}{2(4.94)} t = 0.649

Observe que el resultado tiene sentido, pues según la tabla, el máximo debería estar entre 0.60 y 0.75

b) El cruce por la «barrera», corresponde al desplazamiento del balón en el eje Y a 9 metros: y(t) = 9.
El polinomio modelo de primer grado usa dos puntos para la interpolación, de los sugeridos se escogen dos, por ejemplo: 0.0 y 0.3.

Usando diferencias finitas avanzadas :

d1 = (8.88-0) = 8.88 y(t) = 0 + 8.88\frac{(t-0)}{1!(0.3)} y(t) = 29.6 t

usando y(t) = 9

29.6 t = 9
t = 0.30
z(0.30) = 1.40
(de la tabla de datos)

cuya respuesta con solo dos dígitos decimales es coherente, al estar cercano el valor a una distancia y=8.88 o aproximado a 9.
La respuesta puede ser mejorada usando más digitos significativos en las operaciones.

c)  La desviación máxima en eje X.
Determine un polinomio para la trayectoria en el eje X y obtenga el máximo igualando la derivada del polinomio x(t) a cero.

Por simplicidad, se usa un polinomio de segundo grado, alrededor de los valores máximos en el eje X

t x(t) d1 d2 d3
0.60 1.80 0.20 -0.30
0.75 2.00 -0.10
0.90 1.90
x(t) = 1.80 + 0.20 \frac{(t-0.60)}{1!(0.15)} -0.30 \frac{(t-0.60)(t-0.75)}{2!(0.15)^2} x(t) = 1.80 + 1.33 (t-0.60) - 6.67(t-0.60)(t-0.75)

como se busca el máximo, se usa la derivada \frac{dx}{dt} = 0

\frac{dx}{dt} = 1.33 - 6.67(2t +(-0.60-0.75)) 1.33 - 13.34t + 9.00 = 0 -13.34t + 10.33 = 0

t = 0.77

x(0.77) = 1.80 + 1.33(0.77-0.60) - 6.67(0.77-0.60)(0.77-0.75) x(0.77) = 2.003

lo que es coherente con la tabla para el eje x, pues el máximo registrado es 2, y el próximo valor es menor, la curva será decreciente.


Desarrollo extra para observar y verificar resultados:

Usando los puntos y un graficador 3D se puede observar la trayectoria:

Tarea: Realice el ejercicio usando los algoritmos en Python, muestre los polinomios obtenidos y compare.

Nota: La gráfica 3D presentada usa interpolación de Lagrange con todos los puntos. Realice las observaciones y recomendaciones necesarias y presente su propuesta como tarea.

s1Eva_IT2009_T2_AN Materiales y Productos 3×4

1. Plantear el sistema de ecuaciones:
0.2 x0 + 0.5 x1 + 0.4 x2 + 0.2 x3 = 10
0.3 x0 +  0  x1 + 0.5 x2 + 0.6 x3 = 12
0.4 x0 + 0.5 x1 + 0.1 x2 + 0.2 x3 = 15

Observe que hay más incógnitas que ecuaciones.

Para equiparar las ecuaciones con el número de incógnitas, podriamos suponer que uno de los productos NO se fabricará. por ejempo el producto x3 que podría hacerse igual a cero. Supone que la variable libre es x3 .

Para mantener la forma de las ecuaciones para otros valores de x3, se pasa la variable y su coeficiente a la derecha.

0.2 x0 + 0.5 x1 + 0.4 x2 = 10 - 0.2 x3 
0.3 x0 +  0  x1 + 0.5 x2 = 12 - 0.6 x3 
0.4 x0 + 0.5 x1 + 0.1 x2 = 15 - 0.2 x3 

Para analizar el ejercicio, se supondrá que el valor de x3 = 0, lo que permite usar el modelo del problema como A.X=B .En caso de que x3 sea diferente de cero,  el vector B modifica, y se puede proceder con el sistema de ecuaciones.

2. Convertir a la forma matricial AX = B

Siendo así, suponiendo que x3 = 0, el ejercicio se puede desarrollar usando:

A = np.array([[0.2, 0.5, 0.4],
              [0.3, 0.0, 0.5],
              [0.4, 0.5, 0.1]])
B = np.array([[10],
              [12],
              [15]])

que luego armar el algoritmo y su ejecución, obtendría una solución semejante:

array([[ 32.10526316],
       [  3.36842105],
       [  4.73684211]])

Nota: Para desarrollar, se recomienda seguir los pasos indicados en:

Unidad03 y 04: Sistema de ecuaciones – Pasos

Encontrada la solución, modifique el algoritmo para calcular B en función de x3, pregunte el valor al inicio, y vuelva a calcular.

x3 = float(input('cantidad a producir de cuarto producto: ')

Observe el rango para x3, por ejemplo que debe ser mayor o igual que cero, pues no hay producción negativa. De producir solamente ése un producto, el valor máximo de unidades a obtener no superaría lo posible con la cantidad de materiales disponible.

Ejemplo:

x3 = 1
B3 = np.array([[0.2],
	       [0.6],
	       [0.2]])
Bnuevo = B - x3*B3

>>> Bnuevo
array([[  9.8],
       [ 11.4],
       [ 14.8]])

Y se vuelve a generar el sistema A.X=Bnuevo

Observación: Algunos productos se fabrican por unidades, no necesesariamente por peso o volumen. ¿comentarios al respecto?

 

s1Eva_IIT2011_T2_AN Sistema de Ecuaciones, diagonal dominante

para resolver, se usa:

A = np.array([[-2, 5, 9],
              [ 7, 1, 1],
              [-3, 7,-1]])
B = np.array([1,6,-26])

debiendo notar que el vector B se encuentra en forma de fila, no de columna.
Para crear la matriz aumentada se usa B en forma de columna.

Para el pivoteo parcial por filas, se usa la matriz aumentada.

Para aplicar el algoritmo de Gauss-Seidel desarrollado, se debe actualizar A y B desde la matriz aumentada a la que se aplicó pivoteo parcial.

El resto es semejante al procedimento de Gauss-Seidel, obteniendo la respuesta de:

matriz aumentada con diagonal dominante: 
[[  7   1   1   6]
 [ -3   7  -1 -26]
 [ -2   5   9   1]]
nuevas matrices: 
[[ 7  1  1]
 [-3  7 -1]
 [-2  5  9]]
[  6 -26   1]
respuesta: 
[[ 1.00000245]
 [-2.99999785]
 [ 1.99999935]]
>>> 

Algoritmo desarrollado por partes:

# 1ra Evaluación II Término 2011
# Tema 2. sistema de ecuaciones
# Desarrolla: diagonal estrictamente dominante
#             Gauss-Seidel
# Vector B se ingresa como fila, luego de convierte a columna

import numpy as np

# INGRESO
A = np.array([[-2, 5, 9],
              [ 7, 1, 1],
              [-3, 7,-1]])
B = np.array([1,6,-26])

# PROCEDIMIENTO
# convierte B a columna
B = np.transpose([B])
# matriz aumentada
aumentada = np.concatenate((A,B), axis=1)
tamano = np.shape(aumentada)
n = tamano[0]
m = tamano[1]
# Diagonal estrictamente dominante
i = 0
while not(i>=(n-1)):
    columna = np.abs(aumentada[i:,i])
    dondemax = np.argmax(columna)
    if (dondemax != 0):
        temporal = np.copy(aumentada[i,:])
        aumentada[i,:] = aumentada[dondemax+i,:]
        aumentada[dondemax+i,:] = temporal
    i=i+1

# Separa matriz aumentada para forma estandar A.X=B
ultima = m-1
penultima = m-2
# antes de última
A = aumentada[:,:ultima]
# solo última columna
B = aumentada[:,ultima] 

# Método GAUSS-SEIDEL
# el ejercicio indica vector X inicia con ceros
# puede ser otro vector
X = np.zeros(n, dtype=float)
tolera = 1e-4

# Método Gauss-Seidel
iteramax = 100
tamano = np.shape(A)
n = tamano[0]
m = tamano[1]
diferencia = np.ones(n, dtype=float)
errado = np.max(diferencia)

itera = 0
while (errado>tolera or itera>iteramax):
    for i in range(0,n,1):
        nuevo = B[i]
        for j in range(0,m,1):
            if (i!=j): # excepto diagonal de A
                nuevo = nuevo-A[i,j]*X[j]
        nuevo = nuevo/A[i,i]
        diferencia[i] = np.abs(nuevo-X[i])
        X[i] = nuevo
    errado = np.max(diferencia)
    itera = itera + 1
# Vector en columna
X =  np.transpose([X])
# No converge
if (itera>iteramax):
    X=0

respuesta = X
# SALIDA
print('matriz aumentada con diagonal dominante: ')
print(aumentada)
print('nuevas matrices: ')
print(A)
print(B)
print('respuesta: ')
print(respuesta)





s1Eva_IIT2011_T2_AN-Sistema-de-ecuaciones

Propuesta de solución

matriz aumentada: 
[[ -2   5   9   1]
 [  7   1   1   6]
 [ -3   7  -1 -26]]
diagonal dominante: 
[[  7   1   1   6]
 [ -3   7  -1 -26]
 [ -2   5   9   1]]
nuevas matrices: 
A: 
[[ 7  1  1]
 [-3  7 -1]
 [-2  5  9]]
B: 
[  6 -26   1]
solucion X: 
[ 1.00000245 -2.99999785  1.99999935]
verifica A.X=B
[  6.00001868 -25.99999164   1.        ]
>>> 

Desarrolle manualmente el ejercicio para verificar sus respuestas.
el archivo de matg1013.py con las funciones desarrolladas en clase se encuentra en la seccion de Recursos de estudio.

# 1ra Evaluación II Término 2011
# Tema 2. Sistema de ecuaciones
import numpy as np
import matg1013 as fcnm

A = np.array([[-2, 5, 9],
              [ 7, 1, 1],
              [-3, 7,-1]])
B = np.array([1,6,-26])
tolera = 1e-4
X = np.zeros(len(B), dtype=float)
iteramax = 100

# PROCEDIMIENTO

# Diagonal dominante:
B = np.transpose([B]) # en columna
aumentada = np.concatenate((A,B),axis=1)
dominante = fcnm.pivoteafila(aumentada)
# Separa A y B de dominante
A = dominante[:,:-1]
B = dominante[:,-1]

respuesta = fcnm.gauss_seidel(A,B,X,tolera, iteramax)
verifica = np.dot(A,respuesta)

# SALIDA
print('matriz aumentada: ')
print(aumentada)
print('diagonal dominante: ')
print(dominante)
print('nuevas matrices: ')
print('A: ')
print(A)
print('B: ')
print(B)
print('solucion X: ')
print(respuesta)
print('verifica A.X=B')
print(verifica)