Autor Tema: Método de integracion de Romberg

0 Usuarios y 1 Visitante están viendo este tema.

27 Abril, 2024, 12:24 am
Leído 50 veces

esmeraldabrown

  • $$\Large \color{#5e8d56}\pi\,\pi\,\pi$$
  • Mensajes: 126
  • País: ve
  • Karma: +0/-1
INTEGRACIÓN DE ROMBERG:

He realizado un código para el siguiente ejercicio:

Use la integración de Romberg para calcular las siguientes aproximaciones \( \displaystyle\int_{1}^{48}\sqrt[ ]{1+(cosx)^2}dx \)

1 [Nota: Los resultados de este ejercicio son interesantes si maneja aritmética
entre siete y nueve dígitos].

(a) Determine Ri,1 para\(  i = 1,2,3,4 y 5 \). Utilice estas aproximaciones para predecir el valor de la
integral.

(b) Determine Ri,i para i = \( 1,2,3,4 y 5 \).Modifique su predicción.

(c) Determine Ri,i para \( i = 7,8,9 y 10 \) .Haga una predicción final.

el código es el siguiente:

import numpy as np

# Función a integrar
def f(x):
    return np.sqrt(1 + np.cos(x)**2)

# Regla del trapecio compuesta
def trapezoid_rule(f, a, b, n):
    h = (b - a) / n
    x = np.linspace(a, b, n+1)
    y = f(x)
    return h/2 * (y[0] + 2*np.sum(y[1:n]) + y[n])

# Integración de Romberg
def romberg_integration(f, a, b, n_max):
    R = np.zeros((n_max+1,n_max+1))

    for i in range(1, n_max+1):
        R[1] = trapezoid_rule(f, a, b, 2**(i-1))

        for j in range(2,i+1):
            R[j] = (4**(j-1) * R[j-1] - R[i-1][j-1]) / (4**(j-1) - 1)

    return R

# Calcular aproximaciones Ri,1 para i=1 a 5
R_approximations = romberg_integration(f, 1, 48, 5)

print("Aproximaciones Ri,1:")
for i in range(1,6):
    print("R{},1: {:.9f}".format(i,np.around(R_approximations[1], 9)))

# Predicción del valor de la integral con la última aproximación Ri,i
integral_prediction = R_approximations[5][5]
print("\nPredicción del valor de la integral con la última aproximación Ri,i: {:.9f}".format(np.around(integral_prediction, 9)))

# Calcular aproximaciones Ri,i para i=1 a 5
R_diagonals = [R_approximations for i in range(6)]

print("\nAproximaciones Ri,i:")
for i in range(6):
    print("R{},{}: {:.9f}".format(i,i,np.around(R_diagonals, 9)))

# Calcular aproximaciones Ri,i para i=7 a 10
R_diagonals_extended = [trapezoid_rule(f, 1, 48 ,2**i) for i in range(7,11)]

print("\nAproximaciones extendidas Ri,i:")
for i in range(7,11):
    print("R{},{}: {:.9f}".format(i,i,np.around(R_diagonals_extended[i-7], 9)))

# Predicción final con la última aproximación extendida Ri,i
integral_prediction_final = R_diagonals_extended[-1]
print("\nPredicción final del valor de la integral con la última aproximación extendida Ri,i: {:.9f}".format(np.around(integral_prediction_final, 9))) ,

me pueden decir si estos resultados están correctos? el codigo arroja estos resultados:
Aproximaciones Ri,1:
\( R1,1: 54.613366966//
R2,1: 57.495690027//
R3,1: 57.062006029//
R4,1: 57.058461760//
R5,1: 57.056432956//
 \)
Predicción del valor de la integral con la última aproximación Ri,i: \( 57.055422056 \)

Aproximaciones Ri,i:
\( R0,0: 0.000000000//
R1,1: 54.613366966//
R2,2: 58.456464381//
R3,3: 56.814843384//
R4,4: 57.070598893//
R5,5: 57.055422056 \)

Aproximaciones extendidas Ri,i:
\( R7,7: 57.158871884//
R8,8: 57.158987051//
R9,9: 57.159016703//
R10,10: 57.159024165 \)

Predicción final del valor de la integral con la última aproximación extendida Ri,i: \( 57.159024165 \)