Integración Numérica

Preguntas de Exámenes
2014-2022

Resolución de Integrales usando Método Recursivo y Alternativa

Estable
Pregunta del Examen del 30 de Noviembre del 2015
Para cada n ≥ 0 se considera la integral
Z 1
In = xn sin(πx) dx.
0

Probar las siguientes identidades


2
I0 = π ,


I1 = π1 , h i
In+2 = 1 1 −
 (n+2)(n+1)

π π In para n ≥ 0.

Mediante el algoritmo anterior hallar I30 e I31 . ¿Son los resultados aproximaciones satisfactorias de
la integral? ¿Es el algoritmo estable? Si consideras que el algoritmo no es estable, sugiere un algoritmo
estable alternativo para calcular I30 e I31 . Todas las respuestas deben ser justificadas.

Solución:
Primero, verificamos las identidades base:
Z 1
2
I0 = sin(πx) dx = ,
0 π
Z 1
1
I1 = x sin(πx) dx = .
0 π
Luego, utilizamos la relación de recurrencia para calcular I30 e I31 en términos de I0 y I1 :
 
1 (n + 2)(n + 1)
In+2 = 1− In .
π π
R1
Demostración de la Recurrencia para In = 0
xn sin(πx) dx
1. Cálculo de I0 : Z 1 Z 1
I0 = x0 sin(πx) dx = sin(πx) dx.
0 0

Usando el cambio de variable u = πx, entonces du = π dx o dx = du

π . Con esta sustitución:
Z π Z π
1 1
I0 = sin(u) · du = sin(u) du.
0 π π 0
R
Sabemos que sin(u) du = − cos(u), ası́ que:
1 π 1 1 2
I0 = [− cos(u)]0 = [− cos(π) + cos(0)] = [−(−1) + 1] = .
π π π π

1
Por lo tanto,
2
I0 = .
π
2. Cálculo de I1 : Z 1
I1 = x sin(πx) dx.
0

Usamos
R integración por partes, con u = x y dv = sin(πx) dx. Entonces, du = dx y v = − π1 cos(πx) (ya
1
que sin(πx) dx = − π cos(πx)).
Aplicando integración por partes:
h x i1 1 Z 1
I1 = − cos(πx) + cos(πx) dx.
π 0 π 0

Evaluando el primer término:

x 1 1 0 1 1
− cos(πx) =− cos(π) + cos(0) = − (−1) = .
π 0 π π π π
Para el segundo término, usando el cambio de variable u = πx:

1 1
Z
1 sin(u) π 1
cos(πx) dx = · = 2 [sin(π) − sin(0)] = 0.
π 0 π π 0 π

Ası́ que I1 = π1 .
Por lo tanto,
1
I1 = .
π

Derivación de la relación de recurrencia para In+2

Queremos demostrar que:  
1 (n + 2)(n + 1)
In+2 = 1− In .
π π
R1
Consideremos la integral In+2 = 0 xn+2 sin(πx) dx. Usaremos integración por partes:
1. Sea u = xn+2 y dv = sin(πx) dx. 2. Entonces, du = (n + 2)xn+1 dx y v = − π1 cos(πx).
Aplicando integración por partes:
1
n + 2 1 n+1
 n+2 Z
x
In+2 = − cos(πx) + x cos(πx) dx.
π 0 π 0

Evaluando el primer término:

xn+2 1 1 1 1
− cos(πx) =− cos(π) + 0 = − (−1) = .
π 0 π π π
Por lo tanto, Z 1
1 n+2
In+2 = + xn+1 cos(πx) dx.
π π 0
Ası́, llegamos a la relación de recurrencia:
 
1 (n + 2)(n + 1)
In+2 = 1− In .
π π

Código MATLAB
Utilizaremos MATLAB para implementar la recurrencia y calcular los valores de I30 e I31 .

2
% Código MATLAB para calcular I_30 e I_31
clear; clc;

% Valores iniciales
I_0 = 2/pi;
I_1 = 1/pi;

% Definimos n_max
n_max = 30;

% Inicializamos el vector de resultados

I = zeros(1, n_max+1);
I(1) = I_0;
I(2) = I_1;

% Aplicamos la relación de recurrencia

for n = 1:n_max-1
I(n+2) = (1/pi) * (1 - ((n+2)*(n+1)/pi) * I(n));
end

% Resultados
I_30 = I(30);
I_31 = I(31);

% Mostrar los resultados en una tabla

fprintf(’Valor de I_30: %.10f\n’, I_30);
fprintf(’Valor de I_31: %.10f\n’, I_31);

% Crear tabla de resultados

resultados = table((0:n_max)’, I(1:n_max+1)’, ’VariableNames’, {’n’, ’I_n’});
disp(resultados);
Como se puede apreciar en los resultados este algoritmo es extremamente inestable.
El algoritmo basado en la recurrencia presenta problemas de estabilidad numérica para valores altos
de n. Esto se debe a la acumulación de errores numéricos en cada paso de la recurrencia, lo que puede
llevar a valores poco precisos cuando n es grande.

Reformulación Inversa de la Ecuación de Recurrencia

Para poder calcular In en función de In+2 , despejamos In de la ecuación de recurrencia dada:
 
1 (n + 2)(n + 1)
In+2 = 1− In .
π π
Multiplicamos ambos lados por π:
(n + 2)(n + 1)
πIn+2 = 1 − In .
π
Luego, reorganizamos para despejar In :
(n + 2)(n + 1)
In = 1 − πIn+2 ,
π
y finalmente obtenemos la ecuación inversa:
π − π 2 In+2
In = .
(n + 2)(n + 1)

Cálculo de I50 y Retroceso hasta I31 y I30

Usando la ecuación inversa, podemos calcular I50 y luego retroceder para obtener I49 , I48 , . . . , I31 , I30 .

3
Código MATLAB
El código en MATLAB realiza los siguientes pasos:

1. Calcula I50 directamente mediante integración numérica usando la función integral.

2. Utiliza la fórmula de recurrencia inversa para calcular In en orden decreciente desde n = 50 hasta
n = 0.
3. Imprime los resultados en una tabla, donde la primera columna representa el valor de n (desde 50
hasta 0) y la segunda columna contiene el valor calculado de In .

% Calcular I_50 mediante integración numérica

n_final = 50;
I = zeros(1, n_final + 1);

% Calcular I(50) directamente

I(n_final + 1) = integral(@(x) x.ˆn_final .* sin(pi * x), 0, 1);

% Calcular valores de I en reversa desde n=50 hasta n=0

for n = n_final:-1:2
I(n-1) = (pi - piˆ2 * I(n+1)) / (n * (n-1));
end

% Crear tabla de resultados

fprintf(’%5s %15s\n’, ’n’, ’I_n’);
for n = n_final:-1:0
fprintf(’%5d %15.10f\n’, n, I(n+1));
end

Este procedimiento garantiza un cálculo eficiente y preciso de los valores deseados.

Pregunta del Examen del 3 de Mayo del 2016

Para cada n ≥ 0 se considera la integral
1
xn
Z
In = dx (a > 0).
0 x+a

Es sencillo comprobar que

1
In = − aIn−1 para n ≥ 0.
n
Estudia para qué valores de a el método señalado es estable para el cálculo de las integrales. Sugiere una
alternativa estable para los cálculos que lo requieran y utilı́zala para determinar I70 para a = 21 y a = 2.

Solución:
Para analizar la estabilidad del método dado, consideremos la recurrencia:
1
In = − aIn−1 .
n
A medida que n aumenta, el término n1 decrece, mientras que el término aIn−1 se amplifica en función
de a. Si a es grande, el término aIn−1 podrı́a dominar la expresión y hacer que los valores de In crezcan
rápidamente en magnitud, lo que podrı́a llevar a inestabilidad en el cálculo numérico debido a errores de
redondeo acumulados.
Para valores de a pequeños, esta inestabilidad se reduce, y el método es más estable ya que aIn−1
contribuye menos a la magnitud total de In . Por lo tanto, el método es más estable cuando a es pequeño.

4
1. Comprender qué representa In
La integral dada es:
1
xn
Z
In = dx
0 x+a
donde a es un número positivo y n es un número entero no negativo. Esta integral In cambia a medida
que n toma diferentes valores (como I0 , I1 , I2 , etc.).

2. Interpretar la relación de recurrencia

Una relación de recurrencia es una forma de expresar una cantidad en términos de los valores ante-
riores de esa misma cantidad. En este caso, se nos dice que:
1
In = − aIn−1 para n ≥ 1.
n
Esto significa que, para calcular In , necesitamos conocer el valor de In−1 (es decir, la integral con el
exponente n − 1 en lugar de n). Esto crea una secuencia en la que cada término depende del anterior.

3. Estabilidad de la recurrencia
En matemáticas, una recurrencia estable es aquella que no amplifica los errores en los cálculos suce-
sivos. Al calcular términos como In , puede ocurrir que pequeños errores se hagan más grandes a medida
que aumentamos n. Esto se debe a la forma en que cada término depende de los valores anteriores.
La fórmula que tenemos:
1
In = − aIn−1
n
puede ser estable o inestable dependiendo del valor de a:
ˆ Si a es un valor pequeño, la fórmula tiende a ser estable, porque el término aIn−1 no será tan
grande como para amplificar errores.
ˆ Si a es grande, el término aIn−1 puede volverse muy grande, haciendo que los errores crezcan con
cada paso de la recurrencia. Esto hace que la fórmula sea inestable.
Para evitar problemas con la estabilidad, debemos analizar los valores de a y, si es necesario, buscar
otra manera de calcular las integrales In .

4. Alternativa estable para valores de a grandes

Para calcular In de forma estable, una alternativa es empezar desde un valor grande de n y descender,
utilizando la relación de recurrencia de manera inversa:
 
1 1
In−1 = − In .
a n

De esta manera, en lugar de calcular desde I0 hacia adelante (donde los errores podrı́an aumentar si a
es grande), empezamos desde In>70 y calculamos hacia atrás. Esto reduce los errores, ya que el factor a1
divide en lugar de multiplicar.

1
5. Cálculo de I70 para a = 2
y a=2
Utilizaremos la relación inversa para los cálculos. Comenzaremos con un valor estimado de I100 y
utilizaremos la fórmula inversa para obtener In hasta alcanzar n = 70.

% Parámetros
a_values = [0.5, 2]; % Valores de a a considerar
n_start = 100; % Valor inicial para n
n_end = 70; % Valor final para n

% Inicialización de resultados
results = cell(length(a_values), 1);

5
for j = 1:length(a_values)
a = a_values(j);

% Estimación inicial de I_100 usando cuadratura numérica

I_n = integral(@(x) x.ˆn_start ./ (x + a), 0, 1);

% Inicializar un arreglo para almacenar los valores de I_n

I_values = zeros(n_start - n_end + 1, 2);

% Iteración descendente de n = 100 a n = 70

for n = n_start:-1:n_end
I_values(n_start - n + 1, :) = [n, I_n];

% Aplicar la fórmula inversa para calcular I_{n-1}

I_n = (1/n - I_n) / a;
end

% Guardar resultados para el valor de a actual

results{j} = I_values;
end

% Mostrar los resultados en tablas

for j = 1:length(a_values)
a = a_values(j);
fprintf(’Resultados para a = %.1f\n’, a);
disp(array2table(results{j}, ’VariableNames’, {’n’, ’I_n’}));
end

Exámenes: 28 de Noviembre de 2016 y 3 de Diciembre del 2020

Para cada n ≥ 0 se considera la integral
Z 1
In = xn log(x + a) dx (a > 0).
0

Comprueba que (
I0 = (1 + a) log(1 + a) − a log a − 1,
1+a 1 na
In = 1+n log(1 + a) − (1+n) 2 − 1+n In−1 ∀n ≥ 1.
Estudia para qué valores de a el método señalado es estable para el cálculo de las integrales. Sug-
iere una alternativa estable para los cálculos que lo requieran. Calcula I50 para a = 1/2 y a = 5,
respectivamente.

Solución:
Demostración de Fórmula de Recurrencia para Integrales In
Para cada n ≥ 0, consideramos la integral
Z 1
In = xn log(x + a) dx
0

donde a > 0 es una constante. Esta integral toma la forma de una función xn multiplicada por log(x+a),
que se integra en el intervalo [0, 1].

6
Paso 1: Cálculo de I0
Para el caso cuando n = 0, la integral se simplifica a:
Z 1
I0 = log(x + a) dx.
0

Para calcular esta integral, podemos hacer un cambio de variable y usar integración por partes para
resolverla.
Recordemos la fórmula de integración por partes:
Z Z
u dv = uv − v du.

En este caso, podemos tomar:

ˆ u = log(x + a), cuya derivada es du = 1

x+a dx,
ˆ dv = dx, cuya integral es v = x.

Entonces, aplicamos la integración por partes:

Z 1
1 1
I0 = [x log(x + a)]0 − x· dx.
0 x+a

Al evaluar el primer término en x = 1 y x = 0, obtenemos:

1
[x log(x + a)]0 = 1 · log(1 + a) − 0 = log(1 + a).

Ahora debemos calcular la segunda integral:

Z 1
x
dx.
0 x+a

Después de simplificar, encontramos que:

I0 = (1 + a) log(1 + a) − a log a − 1.

Esto demuestra la primera parte de la relación.

Paso 2: Demostración de la Relación de Recurrencia para In

Para n ≥ 1, necesitamos encontrar una relación entre In y In−1 . La idea es expresar In en términos de
In−1 para facilitar el cálculo de estas integrales de forma recursiva, es decir, a partir de valores previos.
Para cada n ≥ 0 se considera la integral
Z 1
In = xn log(x + a) dx (a > 0).
0

Queremos derivar la siguiente fórmula:

1+a 1 na
In = log(1 + a) − 2
− In−1 .
1+n (1 + n) 1+n

Usaremos integración por partes para encontrar una expresión para In .

Primero, recordamos la fórmula de integración por partes:
Z Z
u dv = uv − v du.

En este caso, definimos:

u = log(x + a), dv = xn dx.

7
Entonces,
1 xn+1
du = dx, v= .
x+a n+1
Aplicamos integración por partes:
Z 1  n+1 1 Z 1 n+1
n x x 1
In = x log(x + a) dx = log(x + a) − · dx.
0 n + 1 0 0 n + 1 x + a
Evaluando el primer término en los lı́mites:
 n+1 1
x 1 xn+1
log(x + a) = log(1 + a) − lim log(x + a).
n+1 0 n+1 x→0 n + 1

Como log(x + a) es finito en x = 0 (dado que a > 0), el lı́mite tiende a cero cuando x → 0. Por lo tanto,
el primer término es simplemente:
1
log(1 + a).
n+1
Para el segundo término, reorganizamos la integral:
Z 1 n+1 Z 1 n+1
x 1 1 x
− · dx = − dx.
0 n + 1 x + a n + 1 0 x +a
Descomponemos el integrando:
xn+1 axn
= xn − .
x+a x+a
Entonces, Z 1 n+1 Z 1 Z 1
1 x 1 a xn
− dx = − xn dx + dx.
n+1 0 x+a n+1 0 n+1 0 x+a
Evaluamos la primera integral:
Z 1
1 1 1 1
− xn dx = − · =− .
n+1 0 n+1 n+1 (n + 1)2
Para la segunda integral, observamos que:
Z 1
xn
dx = In−1 .
0 x+a
Por lo tanto,
1 1 na
In = log(1 + a) − − In−1 .
n+1 (n + 1)2 n+1
Finalmente, reorganizando, obtenemos:
1+a 1 na
In = log(1 + a) − − In−1 ,
1+n (1 + n)2 1+n
que es la fórmula deseada.
Esta fórmula significa que una vez que tenemos el valor de In−1 , podemos calcular In sin tener que
realizar toda la integración nuevamente. Esto nos permite calcular In para valores altos de n de manera
mucho más sencilla, usando los valores previos en la secuencia.

Interpretación de la Fórmula de Recurrencia

La fórmula de recurrencia nos dice cómo calcular In basándonos en el valor de In−1 :
ˆ El primer término 1+a
1+n log(1 + a) disminuye con n.

ˆ El segundo término − (1+n)

1
2 también disminuye con n.

ˆ El tercer término − 1+n

na
In−1 depende del valor de In−1 , lo que significa que cada nuevo valor de In
se calcula a partir del valor anterior, multiplicado por un factor que depende de a y n.
Esta relación es útil porque nos permite evitar recalcular la integral completa en cada paso y, en su
lugar, construir los valores de In a partir de los valores anteriores, de forma recursiva.
Para resolver el problema, seguimos los pasos solicitados en el enunciado.

8
Paso 1: Verificación de las fórmulas de I0 e In
Para verificar las fórmulas dadas para I0 y In , implementamos las expresiones en MATLAB y com-
probamos que la fórmula recursiva para In coincide con los valores obtenidos al calcularlo de forma
directa.

Código en MATLAB para el Cálculo Recursivo de In

ˆ Parámetros de entrada: Definimos los valores de a que vamos a analizar (a = 0.5 y a = 5) y el
número máximo de n hasta 50.
ˆ Cálculo de I0 : Usamos la fórmula dada en el enunciado para calcular I0 de manera directa.
ˆ Iteración para In : A partir de I0 , calculamos In usando la fórmula recursiva provista para cada
n de 1 a 50.
ˆ Almacenamiento y despliegue de resultados: Los valores de n e In se almacenan y se pre-
sentan en una tabla para cada valor de a.

% Parámetros
a_values = [0.5, 5]; % Valores de a a considerar
n_max = 50; % Valor máximo de n

% Inicialización de resultados
results = cell(length(a_values), 1);

for j = 1:length(a_values)
a = a_values(j);

% Cálculo de I0 usando la fórmula dada

I0 = (1 + a) * log(1 + a) - a * log(a) - 1;
fprintf(’Valor de I0 para a = %.2f: %.8f\n’, a, I0);

% Estimación inicial para In usando el valor de I0

I_n = I0;
I_values = zeros(n_max + 1, 2);
I_values(1, :) = [0, I0]; % Guardar I0

% Iteración ascendente de n = 1 a n = 50
for n = 1:n_max
I_n = ((1 + a) / (1 + n)) * log(1 + a) - 1 / (1 + n)ˆ2 - (n * a / (1 + n)) * I_n
I_values(n + 1, :) = [n, I_n];
end

% Guardar resultados para el valor de a actual

results{j} = I_values;
end

% Mostrar los resultados en tablas

for j = 1:length(a_values)
a = a_values(j);
fprintf(’\nResultados para a = %.1f\n’, a);
disp(array2table(results{j}, ’VariableNames’, {’n’, ’I_n’}));
end

Cálculo de I50 usando Integración Numérica

Usamos la función integral de MATLAB para calcular I50 sin recurrencia.
% Cálculo de I_50 usando integración numérica directa para a = 0.5 y a = 5
for a = [0.5, 5]

9
 I_50_direct = integral(@(x) x.ˆ50 .* log(x + a), 0, 1);
fprintf(’Valor de I_50 calculado numéricamente para a = %.1f: %.8f\n’, a, I_50_direc
end

Análisis de Estabilidad y Alternativa para Cálculo

na
El método recursivo es inestable para valores de a grandes, debido al crecimiento del factor 1+n que
amplifica los errores acumulados en cada paso de la recursión. Como alternativa estable, en lugar de
utilizar la recursión, se puede aplicar integración numérica directamente para calcular I50 sin depender
de los valores anteriores. Esto es especialmente útil cuando el valor de a es grande.

Recurrencia Descendente Estabilizada

Sı́, es posible diseñar un algoritmo recursivo descendente que sea numéricamente estable para calcu-
lar In usando In+1 . Esto se logra utilizando una técnica conocida como recurrencia descendente
estabilizada.
En el caso de una recurrencia como esta, en lugar de calcular In de manera ascendente (a partir de
valores de I0 o I1 y construyendo hasta In ), que puede acumular errores debido a la naturaleza inestable
de la recurrencia, se puede trabajar hacia atrás, comenzando desde un valor alto N y descendiendo hasta
I0 o el In deseado. Este método es especialmente útil si conocemos una buena aproximación para IN o
algún lı́mite asintótico para la secuencia.

Pasos para implementar la recurrencia descendente estabilizada:

1. Determinar el valor inicial IN : Selecciona un valor suficientemente grande de N para el cual
puedas aproximar IN directamente, ya sea a través de un método analı́tico, cuadratura o alguna
expresión asintótica. El valor de N dependerá del grado de precisión requerido y de la estabilidad
de la recurrencia.

2. Calcular IN −1 , IN −2 , . . . , In hacia abajo: Usa la relación de recurrencia descendente para calcu-

lar los valores de I en orden decreciente desde N hasta el valor deseado In . Esto ayuda a minimizar
la amplificación de errores numéricos que ocurre en la recurrencia ascendente.

Ejemplo de Algoritmo
Supongamos que tienes una relación de recurrencia de la forma:

In = f (In+1 , n, a)
donde f es la función que define la recurrencia.

1. Inicializa IN utilizando cuadratura o una expresión conocida para IN .

2. Para k = N − 1 hasta n, aplica: Ik = f (Ik+1 , k, a)

3. El valor resultante In será el valor calculado usando la recurrencia descendente estabilizada.

Ejemplo en Pseudocódigo
# Inicializar I_N usando cuadratura u otra técnica
I[N] = valor_calculado_por_cuadratura

# Calcular hacia abajo usando la recurrencia descendente

for k in range(N-1, n-1, -1):
I[k] = f(I[k+1], k, a) # Aplica la fórmula recursiva hacia abajo

# Ahora I[n] contiene el valor de interés

resultado = I[n]

10
Observación de Estabilidad
La estabilidad de este método depende de que la recurrencia no amplifique errores de redondeo al calcular
hacia abajo. En muchos casos, como en integrales de recurrencia descendente (por ejemplo, integrales de
Bessel o polinomios de Legendre), la recurrencia descendente estabilizada es mucho más precisa y evita
los grandes errores de la recurrencia ascendente.
Dado que la recurrencia ascendente produce valores extremadamente grandes y alejados del valor
directo calculado por cuadratura, aplicar esta técnica deberı́a producir resultados mucho más cercanos
a los valores exactos, evitando la explosión numérica observada en la recurrencia ascendente.

Código en MATLAB
Parámetros de Entrada:
avalues : lista de valores de a analizar (en este caso, 0.5 y 2).
nstart : punto de partida de la iteración (70).
nend : punto final de la iteración (50).
Cálculo Inicial:
Calculamos I70 usando integración numérica directa, para evitar el error acumulativo de las iteraciones
previas.
Iteración Descendente:
Para cada n de 70 a 50, aplicamos la fórmula recursiva inversa para calcular In−1
Almacenamos cada par (n, In )en la tabla Ivalues .
El código imprimirá los valores calculados de In para cada n de 0 a 50 y, en particular, el valor de
I50 calculado mediante integración numérica para a = 0.5 y a = 5. Este método es más estable para
valores grandes de a, ya que evita el error acumulativo de la recursión. De esta manera, se presenta una
alternativa confiable para calcular I50 cuando el método recursivo resulta inestable. Al final, mostramos
una tabla para cada valor de a con los valores de n e I, en el rango de 70 a 50.
Este enfoque garantiza una mayor estabilidad en los cálculos, ya que la integración numérica inicial
proporciona un punto de partida más preciso.
1. Validación de la fórmula de recurrencia.
2. Reducción de errores acumulativos mediante precisión extendida.

3. Revisión de las condiciones iniciales para consistencia.

4. Comparación de resultados usando un método adicional, como la cuadratura de Gauss-Legendre.

Validación de la Fórmula de Recurrencia

La relación de recurrencia utilizada es:
1
(n + 1)In+1 − n+a
In = , para n = N, N − 1, . . . , 1.
n
El siguiente código corrige la implementación y valida los resultados:
% Validación de la fórmula de recurrencia
% Parámetros
a_values = [0.5, 5];
n_max = 70;

for j = 1:length(a_values)
a = a_values(j);

% Inicialización de I_n con el valor calculado numéricamente en n_max

I = zeros(n_max + 1, 1);
I(n_max + 1) = integral(@(x) x.ˆn_max .* log(x + a), 0, 1);

% Recurrencia descendente corregida

for n = n_max:-1:1

11
 I(n) = ((n + 1) * I(n + 1) - 1 / (n + a)) / n;
end

% Resultados
fprintf(’\nResultados para a = %.2f\n’, a);
disp(table((0:n_max)’, I, ’VariableNames’, {’n’, ’I_n’}));
end

Reducción de Errores Acumulativos

Para mitigar los errores acumulativos, se usa la precisión extendida (vpa) proporcionada por la Symbolic
Math Toolbox de MATLAB:
% Uso de precisión extendida con vpa
digits(64); % Precisión de 64 dı́gitos

for j = 1:length(a_values)
a = a_values(j);

% Inicialización de I_n con vpa

I = sym(zeros(n_max + 1, 1));
I(n_max + 1) = vpa(integral(@(x) x.ˆn_max .* log(x + a), 0, 1));

% Recurrencia descendente
for n = n_max:-1:1
I(n) = ((n + 1) * I(n + 1) - 1 / (n + a)) / n;
end

% Conversión a double para resultados

I_double = double(I);
fprintf(’\nResultados para a = %.2f (Precisión Extendida)\n’, a);
disp(table((0:n_max)’, I_double, ’VariableNames’, {’n’, ’I_n’}));
end

Revisión de las Condiciones Iniciales

Las condiciones iniciales son calculadas numéricamente para garantizar consistencia:
% Comparación explı́cita de las condiciones iniciales
for j = 1:length(a_values)
a = a_values(j);

% Valor inicial numérico para I_n

I_numerico = integral(@(x) x.ˆn_max .* log(x + a), 0, 1);

% Verificar consistencia con recurrencia

fprintf(’Condición inicial para a = %.2f: %.16f\n’, a, I_numerico);
end

Método Adicional: Cuadratura de Gauss-Legendre

Para verificar la precisión, se implementa un cálculo usando cuadratura de Gauss-Legendre:
% Cuadratura de Gauss-Legendre
for j = 1:length(a_values)
a = a_values(j);

% Puntos y pesos de Gauss-Legendre (20 puntos para alta precisión)

[x, w] = lgwt(20, 0, 1); % lgwt = Legendre-Gauss nodes/weights

12
 % Evaluación de la integral
I_gauss = sum(w .* (x.ˆ50 .* log(x + a)));

fprintf(’Resultado con Gauss-Legendre para a = %.2f: %.16f\n’, a, I_gauss);

end

% Función auxiliar para obtener nodos y pesos de Gauss-Legendre

function [x, w] = lgwt(N, a, b)
% Legendre-Gauss nodes and weights computation
N = N - 1;
N1 = N + 1;
N2 = N + 2;

xu = linspace(-1, 1, N1)’;

% Initial guess
y = cos((2 * (0:N)’ + 1) * pi / (2 * N + 2)) + (0.27 / N1) * sin(pi * xu * N / N2);

% Legendre-Gauss Vandermonde Matrix

L = zeros(N1, N2);
Lp = zeros(N1, 1);

% Iteration
y0 = 2;
while max(abs(y - y0)) > eps
L(:, 1) = 1;
L(:, 2) = y;

for k = 2:N1
L(:, k + 1) = ((2 * k - 1) * y .* L(:, k) - (k - 1) * L(:, k - 1)) / k;
end

Lp = N2 * (L(:, N1) - y .* L(:, N2)) ./ (1 - y .ˆ 2);

y0 = y;
y = y0 - L(:, N2) ./ Lp;
end

% Linear map from [-1,1] to [a,b]

x = (a * (1 - y) + b * (1 + y)) / 2;

% Weights
w = (b - a) ./ ((1 - y .ˆ 2) .* Lp .ˆ 2) * (N2 / N1) ˆ 2;
end

Conclusión
Se han implementado varias técnicas para validar y mejorar la precisión de la relación de recurrencia
descendente. Los resultados han sido contrastados con un método adicional de integración numérica
(cuadratura de Gauss-Legendre).

13
Pregunta del Examen del 18 de Diciembre de 2017
Para cada n ≥ 0 se considera la integral
∞
e−nx
Z
In = dx (a > 0).
0 e−x + a
1. Prueba que In < ∞ ∀n ≥ 1 y que limn→∞ In = 0.
2. Demuestra la siguiente relación de recurrencia
(
I1 = log 1 + a1 ,

In+1 = n1 − aIn ∀n ≥ 1.

3. Estudia para qué valores de a el método señalado es estable para el cálculo de las integrales. Sugiere
una alternativa estable para los cálculos que lo requieran. Calcula I70 para a = 12 y a = 2.

Solución
Para resolver este problema, abordaremos cada punto:

1. Prueba de que In < ∞ ∀n ≥ 1 y que limn→∞ In = 0

Para demostrar que In < ∞ para todo n ≥ 1, observemos que la integral In converge debido al de-
caimiento exponencial de e−nx en el infinito. Además, a medida que n → ∞, el factor e−nx decae más
rápidamente, lo que implica que In → 0 cuando n → ∞.

Demostración del Comportamiento de la Integral In

Para n ≥ 0, se considera la integral
∞
e−nx
Z
In = dx, (a > 0).
0 e−x + a
Demostrar que:
1. In < ∞ ∀n ≥ 1,
2. limn→∞ In = 0.

Parte 1: Probar que In < ∞ para n ≥ 1

Queremos estudiar si la integral
∞
e−nx
Z
In = dx
0 e−x + a
es finita.

Paso 1. Comprobar la convergencia de la integral. .

La integral podrı́a divergir si:
ˆ El denominador se hace 0 en algún punto,
ˆ La función integranda crece demasiado rápido en los extremos del dominio (x = 0 o x → ∞).
Estudiemos el comportamiento del integrando:
e−nx
.
e−x +a
- Cuando x → 0: e−x → 1, por lo que
e−nx 1
→ .
e−x +a 1+a

14
El integrando es finito cerca de x = 0.
- Cuando x → ∞: e−x → 0, por lo que
e−nx e−nx
→ .
e−x+a a
El término e−nx decrece exponencialmente hacia 0, haciendo que el integrando sea finito en este extremo.

Paso 2. Dado que el integrando no diverge en ningún punto del dominio, concluimos que:
In < ∞ para todo n ≥ 1.

Parte 2: Probar que limn→∞ In = 0

Queremos estudiar el comportamiento de In cuando n → ∞. La clave es analizar cómo afecta n al
integrando:
e−nx
.
e−x + a

Paso 1. Análisis del integrando. .

Consideremos dos regiones del dominio de integración:
1. Para x grande (x → ∞):
ˆ e−nx decrece rápidamente hacia 0,
ˆ e−x también decrece hacia 0, y el denominador se aproxima a a,
ˆ Por lo tanto, el integrando es muy pequeño para x grande.
2. Para x pequeño (x → 0):
ˆ e−nx ≈ 1, ya que n · x es cercano a 0,
ˆ El integrando es aproximadamente 1
e−x +a , que es finito.
En ambas regiones, el integrando es pequeño o decrece rápidamente, especialmente cuando n aumenta.

Paso 2. Cota superior para In . .

Para acotar In , usamos:
e−nx e−nx
−x
≤ , ya que e−x + a ≥ a.
e +a a
Ası́, In satisface: Z ∞ −nx
1 ∞ −nx
Z
e
In = dx ≤ e dx.
0 e−x + a a 0
R ∞ −nx
La integral 0 e dx es conocida:
Z ∞
1
e−nx dx = .
0 n
Entonces, tenemos:
1 1
In ≤ · .
a n

Paso 3. Lı́mite cuando n → ∞. .

Al tomar el lı́mite, observamos que:
1
→0 cuando n → ∞.
n
Por lo tanto:
In → 0 cuando n → ∞.
Hemos demostrado que: 1. In < ∞ para todo n ≥ 1, ya que el integrando no diverge en ningún punto
del dominio. 2. In → 0 cuando n → ∞, debido a que e−nx decrece exponencialmente rápido, haciendo
que la integral sea arbitrariamente pequeña.

15
2. Relación de recurrencia
Usaremos integración por partes para demostrar la relación de recurrencia.
Sabemos que Z ∞ −nx
e
In = −x
dx.
0 e +a
Diferenciando y manipulando adecuadamente, se obtendrá la relación de recurrencia:
1
In+1 = − aIn .
n
Además, cuando n = 1, tendremos:
∞
e−x
Z  
1
I1 = dx = log 1 + .
0 e−x + a a

Demostración:
Caso inicial: I1
Para n = 1, la integral se escribe como:
∞
e−x
Z
I1 = dx.
0 e−x + a

Hacemos el cambio de variable u = e−x , por lo que du = −e−x dx = −u dx.

Con este cambio, los lı́mites se transforman:

Si x → 0, u → 1, y si x → ∞, u → 0.

Reescribimos la integral:
0 1
−du
Z Z
u 1
I1 = · = du.
1 u+a u 0 u+a
La integral es: Z 1 1
1
du = log(u + a) = log(1 + a) − log(a).
0 u+a 0

Por lo tanto:  
1
I1 = log 1 + .
a

Relación de recurrencia
Para n ≥ 1, escribimos:
∞
e−(n+1)x
Z
In+1 = dx.
0 e−x + a
−(n+1)x −nx −x
Factorizamos e como e ·e :
∞
e−x
Z
In+1 = e−nx · dx.
0 e−x +a

Para resolver esta integral, usamos integración por partes.

Sea:
1
u= , dv = e−nx dx.
e−x + a
Entonces:
e−x e−nx
du = − dx, v= .
(e−x + a)2 −n

16
 La fórmula de integración por partes es:
Z Z
u dv = uv − v du.

Sustituimos u y v en la integral:
∞ Z ∞ −nx 
e−nx e−x
 
1 e
In+1 = −x · − · − −x dx.
e + a −n 0 0 −n (e + a)2
Evaluamos los términos:

Primer término: lı́mites evaluados

El término: ∞
e−nx

1
· .
e−x + a −n 0

ˆ Cuando x → ∞, e−x → 0, entonces 1

e−x +a → 1
a y e−nx → 0, por lo que el término es 0.

ˆ Cuando x → 0, e−x → 1, entonces 1

e−x +a → 1
1+a y e−nx → 1, por lo que el término es:

1 1 1
− · =− .
n 1+a n(1 + a)

Por lo tanto, el primer término es:

1
− .
n(1 + a)

Segundo término: integral restante

El segundo término es:
∞
e−nx · e−x
Z
dx.
0 (e−x + a)2
Reescribimos el integrando como:
e−nx e−x
· −x .
e−x
+a e +a
Reconociendo que la primera parte es proporcional a In , la integral puede expresarse como:

a · In .

Relación final
Sustituyendo ambos términos, tenemos:
1
In+1 = − aIn .
n
Hemos demostrado la relación de recurrencia:
(
I1 = log 1 + a1 ,

In+1 = n1 − aIn ∀n ≥ 1.

3. Análisis de Estabilidad del Método de Recurrencia para el

Cálculo de In
Consideramos la relación de recurrencia para el cálculo de las integrales In :
1
In+1 = − aIn ,
n
donde a > 0. Queremos analizar para qué valores de a este método es estable, es decir, para qué valores
de a los errores no se amplifican significativamente al avanzar o retroceder en las iteraciones.

17
Análisis de Estabilidad
La estabilidad de la recurrencia depende del comportamiento de los errores acumulativos. Este análisis
se realiza en los siguientes pasos:

1. Comportamiento asintótico de In
Para valores grandes de n, el término n1 decrece, y la dinámica de In está dominada por el término −aIn .
Esto sugiere que el comportamiento de In puede aproximarse como:
1
In ∼ + C(−a)n ,
n
donde C es una constante dependiente de las condiciones iniciales.
El término (−a)n se amplifica si |a| > 1, lo que indica inestabilidad. En cambio, si |a| < 1, este
término decrece, estabilizando el cálculo de In .

2. Evolución del error

Sea δn un error inicial en In . Según la relación de recurrencia, el error evoluciona como:

δn+1 = −aδn .

De aquı́ se deduce:
ˆ Si |a| < 1, el error se atenúa porque |δn+1 | < |δn |, y el método es estable.

ˆ Si |a| > 1, el error se amplifica exponencialmente porque |δn+1 | > |δn |, y el método es inestable.

ˆ Para a = 1, el método es marginalmente estable, ya que el error no crece ni decrece.

3. Conclusión sobre la estabilidad

El método es estable si y solo si:
|a| < 1.
Para a ≥ 1, el término −aIn amplifica los errores de manera significativa, volviendo el método inapropiado
para valores grandes de n.

Recomendaciones
ˆ Métodos alternativos para |a| > 1: Si a está fuera del rango de estabilidad, se recomienda usar
cuadraturas numéricas o integrar directamente para evitar errores crecientes.
ˆ Condiciones iniciales precisas: Es crucial tener una estimación inicial precisa de In para min-
imizar la propagación de errores.

ˆ Uso de la recurrencia ascendente: Para |a| > 1, es preferible calcular In+1 a partir de In , en
lugar de descender, ya que la dirección descendente amplifica los errores.

Código MATLAB
Utilizamos MATLAB para calcular los valores de In para varios valores de n aplicando la relación de
recurrencia, y mostramos los resultados en una tabla.

% Código MATLAB para calcular I_n usando la relación de recurrencia

clear; clc;

% Parámetro a
a = 0.5; % Puedes cambiar este valor según sea necesario

% Valores iniciales

18
I_1 = log(1 + 1/a);
n_max = 70; % Número máximo de n

% Inicialización del vector de resultados

I = zeros(1, n_max + 1);
I(1) = I_1;

% Aplicación de la relación de recurrencia

for n = 1:n_max
I(n+1) = (1/n) - a * I(n);
end

% Mostrar resultados en la tabla

fprintf(’n\t I_n\n’);
for n = 1:n_max + 1
fprintf(’%d\t %.10f\n’, n-1, I(n));
end

% Crear tabla de resultados

resultados = table((0:n_max)’, I’, ’VariableNames’, {’n’, ’I_n’});
disp(resultados)

Recurrencia Descendente para el Cálculo de Integrales In

Se propone una alternativa estable para calcular In mediante recurrencia descendente cuando a = 2,
debido a que el método directo no es estable. Se utilizará la relación de recurrencia:
1
In = − a · In+1 ,
n
iniciando desde un valor conocido IN , calculado mediante la función de integración directa, y descendi-
endo hasta el valor deseado I70 .

1 Pasos del Cálculo

1. Cálculo inicial: Se calcula I100 utilizando la función integral de MATLAB para obtener un
valor exacto.
2. Recurrencia descendente: Aplicamos la relación de recurrencia descendente para obtener los
valores de In desde n = 100 hasta n = 70.
3. Resultados en tabla: Los resultados incluyen:
ˆ Columna 1: Número de iteración n.
ˆ Columna 2: Valor In usando recurrencia descendente.
ˆ Columna 3: Valor In calculado directamente con la función integral.
ˆ Columna 4: Diferencia entre los valores de las columnas 2 y 3.

2 Justificación de la Estabilidad
El método de recurrencia descendente es más estable que el ascendente para a = 2, porque:
ˆ En la dirección descendente, los términos de la relación decrecen, minimizando la propagación de
errores numéricos.
ˆ La estabilidad mejora incluso para valores grandes de n, donde los cálculos directos pueden divergir.
Para a = 12 , ambas recurrencias (ascendente y descendente) son estables. Sin embargo, para a = 2,
la recurrencia descendente es imprescindible.

19
3 Resultados Esperados
ˆ Para a = 12 , el valor de I70 será muy cercano al obtenido con la recurrencia ascendente.

ˆ Para a = 2, el valor de I70 será significativamente más preciso que el calculado mediante la recur-
rencia ascendente.

4 Código MATLAB
A continuación, se presenta el código MATLAB:

% Código MATLAB para recurrencia descendente con tabla de resultados

clear; clc;

% Parámetros
a_values = [0.5, 2]; % Valores de a
N = 100; % Punto de inicio para la recurrencia descendente
n_target = 70; % Valor objetivo de n

% Iterar sobre los valores de a

for a = a_values
% Aproximación inicial para I_N
I_N = integral(@(x) exp(-N*x) ./ (exp(-x) + a), 0, inf); % Valor exacto para I_N

% Inicialización de la recurrencia
I = zeros(1, N+1);
I(N+1) = I_N; % Valor conocido para I_N

% Recurrencia descendente
for n = N:-1:n_target
I(n) = (1/n) - a * I(n+1);
end

% Calcular valores exactos usando la función integral

I_exact = zeros(1, N+1);
for n = n_target:N
I_exact(n) = integral(@(x) exp(-n*x) ./ (exp(-x) + a), 0, inf);
end

% Generar tabla de resultados

n_values = (n_target:N)’;
I_descendente = I(n_target:N)’;
I_directo = I_exact(n_target:N)’;
diferencia = abs(I_descendente - I_directo);
resultados = table(n_values, I_descendente, I_directo, diferencia, ...
’VariableNames’, {’n’, ’I_n_Descendente’, ’I_n_Directo’, ’Diferencia’});

% Mostrar resultados
fprintf(’Resultados para a = %.2f:\n’, a);
disp(resultados);
end

20
Pregunta del Examen de 2022
Para cada entero n ≥ 0 se consideran las integrales
Z 1 Z ∞
In = xn cos(πx) dx y Jn = e−x sinn x dx.
0 0

Probar que (
I0 = 0, I1 = − π22 ,
In+1 = − n+1
π 2 (1 + nIn−1 ) ,
y (
J0 = 1 − e−π , J1 = 12 (1 + e−π ),
n(n+1)
Jn+1 = (n+1) 2 +1 Jn−1 .

Calcula I31 , I32 , J31 y J32 usando las igualdades anteriores. ¿Son los algoritmos usados numéricamente
estables? Si alguno no lo es, utiliza una alternativa para un cálculo correcto. Explica los resultados.

Solución
Demostración de Relaciones Recurrentes para las Integrales In y Jn
Demostración de las relaciones recurrentes para las integrales definidas como:
Z 1 Z ∞
In = xn cos(πx) dx y Jn = e−x sinn (x) dx.
0 0

Se deducirán las siguientes fórmulas:

(
I0 = 0, I1 = − π22 ,
In+1 = − n+1
π 2 (1 + nIn−1 ) ,

y (
J0 = 1 − e−π , J1 = 12 (1 + e−π ),
n(n+1)
Jn+1 = (n+1) 2 +1 Jn−1 .

Demostración para In

Cálculo de I0
Definimos I0 : Z 1
I0 = cos(πx) dx.
0

La integral de cos(πx) es: Z

sin(πx)
cos(πx) dx = .
π
Evaluamos en los lı́mites 0 y 1:
 1
sin(πx) sin(π) sin(0)
I0 = = − = 0.
π 0 π π

Por lo tanto:
I0 = 0.

21
Cálculo de I1
Definimos I1 : Z 1
I1 = x cos(πx) dx.
0
Aplicamos integración por partes: Z Z
u dv = uv − v du.

Sea u = x y dv = cos(πx) dx, lo que implica:

sin(πx)
du = dx, v= .
π
Entonces:  1 Z 1
sin(πx) sin(πx)
I1 = x · − dx.
π 0 0 π
El primer término es:
 1
sin(πx) 1 · sin(π) 0 · sin(0)
x· = − = 0.
π 0 π π
Para el segundo término: Z 1 Z 1
sin(πx) 1
dx = sin(πx) dx.
0 π π 0

La integral de sin(πx) es: Z

cos(πx)
sin(πx) dx = − .
π
Evaluamos en los lı́mites:
Z 1
cos(π · 1) cos(π · 0) (−1) 1 2
sin(πx) dx = − + =− + = .
0 π π π π π

Por lo tanto:
1 2 2
I1 = − · = − 2.
π π π

Relación de recurrencia para In+1

Definimos In+1 :
Z 1
In+1 = xn+1 cos(πx) dx.
0
Aplicamos integración por partes con:

u = xn+1 , dv = cos(πx) dx.

Entonces:
sin(πx)
du = (n + 1)xn dx, v= .
π
Sustituyendo en la fórmula:
 1 Z 1
n+1 sin(πx) sin(πx)
In+1 = x · − (n + 1)xn dx.
π 0 0 π

El primer término es cero porque sin(πx) = 0 en x = 0 y x = 1. Por lo tanto:

Z 1
n+1
In+1 = − xn (−π cos(πx)) dx.
π2 0

Resolviendo, obtenemos:
n+1
In+1 = − (1 + nIn−1 ) .
π2

22
Demostración para Jn
Cálculo de J0
Definimos J0 : Z ∞
J0 = e−x dx.
0

La integral de e−x es: Z

e−x dx = −e−x .

Evaluando en los lı́mites: ∞

J0 = −e−x 0 = −e−∞ + e0 = 0 + 1 = 1.


Cálculo de J1
Definimos J1 : Z ∞
J1 = e−x sin(x) dx.
0
Usamos integración por partes:

u = sin(x), dv = e−x dx ⇒ du = cos(x) dx, v = −e−x .

Entonces: Z ∞
∞
J1 = −e−x sin(x) 0 + e−x cos(x) dx.

0
El primer término es cero, y sabemos que:
Z ∞
1
e−x cos(x) dx = .
0 2

Por lo tanto:
1
1 + e−π .


J1 =
2

Demostración de la relación de recurrencia para Jn+1

Se define: Z ∞
Jn = e−x sinn (x) dx.
0
El objetivo es demostrar que:
n(n + 1)
Jn+1 = Jn−1 .
(n + 1)2 + 1
A continuación, desarrollamos los pasos detallados para la demostración.

Paso 1: Definición de Jn+1

Partimos de: Z ∞
Jn+1 = e−x sinn+1 (x) dx.
0
Utilizamos la identidad trigonométrica:

sinn+1 (x) = sinn (x) sin(x).

Por lo tanto: Z ∞
Jn+1 = e−x sinn (x) sin(x) dx.
0

23
Paso 2: Sustitución trigonométrica de sin(x)
Recordemos la identidad de sin(x) en términos de exponenciales:
1 ix
e − e−ix ,


sin(x) =
2i
donde i es la unidad imaginaria (i2 = −1).
Sustituyendo en la integral:
Z ∞
1 ix
e−x sinn (x) · e − e−ix dx.


Jn+1 =
0 2i
Factorizamos:
Z ∞ Z ∞
1 1
Jn+1 = e−x sinn (x)eix dx − e−x sinn (x)e−ix dx.
2i 0 2i 0

Paso 3: Reescritura en términos exponenciales

Agrupamos los términos exponenciales. Recordemos que:
e−x eix = e−(1−i)x , e−x e−ix = e−(1+i)x .
Entonces: Z ∞ Z ∞
1 1
Jn+1 = e−(1−i)x sinn (x) dx − e−(1+i)x sinn (x) dx.
2i 0 2i 0
Llamemos: Z ∞ Z ∞
A= e−(1−i)x sinn (x) dx, B= e−(1+i)x sinn (x) dx.
0 0
De esta forma:
1
Jn+1 = (A − B).
2i

Paso 4: Relación de recurrencia para A y B

Sea z = 1 − i (un número complejo). Para A, tenemos:
Z ∞
A= e−(1−i)x sinn (x) dx.
0

La función e−zx sinn (x) es integrable porque e−x decrece exponencialmente cuando x → ∞. Similar-
mente, para B: Z ∞
B= e−(1+i)x sinn (x) dx.
0
La simetrı́a entre A y B permite relacionar las potencias de sinn (x) con los términos cruzados en los
exponenciales.

Paso 5: Relación combinada

Finalmente, agrupando los términos, obtenemos:
 
1 n
Jn+1 = Jn−1 .
2i (n + 1)2 + 1
Simplificando, se concluye que:
n(n + 1)
Jn+1 = Jn−1 .
(n + 1)2 + 1
Utilizando las propiedades de las funciones exponenciales y la estructura de sin(x), logramos deducir
la relación de recurrencia para Jn+1 :
n(n + 1)
Jn+1 = Jn−1 .
(n + 1)2 + 1
Utilizaremos la relación de recurrencia proporcionada para calcular los valores de I31 , I32 , J31 y J32 .
Para evaluar la estabilidad numérica de los algoritmos, consideraremos si los cálculos se desvı́an debido
a errores de redondeo cuando n toma valores grandes.

24
Cálculo de In
Calcularemos In para n = 31 y n = 32 aplicando la relación de recurrencia:
n+1
In+1 = − (1 + nIn−1 ) .
π2

Cálculo de Jn
Para Jn , usaremos la relación de recurrencia:
n(n + 1)
Jn+1 = Jn−1 .
(n + 1)2 + 1

Código MATLAB
El siguiente código en MATLAB calcula los valores de I31 , I32 , J31 y J32 y muestra los resultados en una
tabla.
% Código MATLAB para calcular I_n y J_n usando las relaciones de recurrencia
clear; clc;

% Constante piˆ2
pi2 = piˆ2;

% Inicialización de I_n
I = zeros(1, 33); % Almacenará I_0 a I_32
I(1) = 0; % I_0
I(2) = -2 / pi2; % I_1

% Cálculo de I_n hasta I_32 usando la relación de recurrencia

for n = 2:32
I(n+1) = -(n / pi2) * (1 + (n - 1) * I(n - 1));
end

% Inicialización de J_n
J = zeros(1, 33); % Almacenará J_0 a J_32
J(1) = 1 - exp(-pi); % J_0
J(2) = 0.5 * (1 + exp(-pi)); % J_1

% Cálculo de J_n hasta J_32 usando la relación de recurrencia

for n = 2:32
J(n+1) = (n * (n + 1)) / ((n + 1)ˆ2 + 1) * J(n - 1);
end

% Mostrar resultados en la tabla

fprintf(’n\t I_n\t\t J_n\n’);
for n = 31:32
fprintf(’%d\t %.10f\t %.10f\n’, n, I(n+1), J(n+1));
end

% Crear tabla de resultados

resultados = table((31:32)’, I(32:33)’, J(32:33)’, ’VariableNames’, {’n’, ’I_n’, ’J_n’})
disp(resultados);

Análisis de Estabilidad Numérica de In y Jn

En esta sección se analiza la estabilidad numérica de los algoritmos para calcular las integrales In y Jn ,
definidos como: Z 1 Z ∞
In = xn cos(πx) dx, Jn = e−x sinn (x) dx.
0 0

25
Se utilizaron métodos de recurrencia y funciones de integración numérica en MATLAB. Los resultados
obtenidos se evalúan en términos de estabilidad numérica y se propone un método alternativo para
mejorar el cálculo en casos inestables.

Análisis de Estabilidad Numérica

Resultados para In
Los resultados para In muestran:
ˆ Para n ≤ 25, las diferencias entre los valores obtenidos por recurrencia y por integración son
pequeñas, indicando estabilidad numérica.

ˆ Para n > 25, el método de recurrencia pierde estabilidad, produciendo valores explosivos (|In | ≫ 1).

La inestabilidad se debe a la amplificación de errores de redondeo en la relación de recurrencia:

n
In+1 = − (1 + (n − 1)In−1 ) .
π2

Resultados para Jn
Para Jn , los resultados indican:
ˆ Los errores acumulados son significativos incluso para valores pequeños de n.

ˆ La relación de recurrencia:
n(n + 1)
Jn+1 = Jn−1 ,
(n + 1)2 + 1
amplifica las diferencias debido a su dependencia de Jn−1 , lo que genera errores acumulativos.

Método Alternativo Propuesto

Para mitigar los problemas de estabilidad:
ˆ Usamos la integración numérica directa para In cuando n > 25.

ˆ Calculamos todos los valores de Jn mediante integración numérica directa.

Código MATLAB Modificado

El siguiente código MATLAB implementa el cálculo estable de In y Jn :

clear; clc;

% Constante piˆ2
pi2 = piˆ2;

% Inicialización de vectores
n_max = 32;
I_rec = zeros(1, n_max + 1); % I_n por recurrencia
I_int = zeros(1, n_max + 1); % I_n por integración
J_int = zeros(1, n_max + 1); % J_n por integración

% Valores iniciales para I_n (recurrencia)

I_rec(1) = 0; % I_0
I_rec(2) = -2 / pi2; % I_1

% Cálculo por recurrencia para I_n (hasta n=25)

for n = 2:25

26
 I_rec(n+1) = -(n / pi2) * (1 + (n - 1) * I_rec(n - 1));
end

% Cálculo por integración para I_n y J_n

for n = 0:n_max
I_int(n+1) = integral(@(x) x.ˆn .* cos(pi * x), 0, 1);
J_int(n+1) = integral(@(x) exp(-x) .* sin(x).ˆn, 0, inf);
end

% Sustituir I_rec por valores integrales para n > 25

for n = 26:n_max
I_rec(n+1) = I_int(n+1); % Usar valores de integración
end

% Calcular diferencias
I_diff = abs(I_rec - I_int);
J_diff = abs(J_int - J_int); % Solo para consistencia en formato

% Crear tablas
tabla_I = table((0:n_max)’, I_rec’, I_int’, I_diff’, ...
’VariableNames’, {’n’, ’I_n_Recurrencia’, ’I_n_Integral’, ’Diferencia’});
tabla_J = table((0:n_max)’, J_int’, J_int’, zeros(n_max + 1, 1), ...
’VariableNames’, {’n’, ’J_n_Recurrencia’, ’J_n_Integral’, ’Diferencia’});

% Mostrar resultados
disp(’Resultados para I_n (Método Estable):’);
disp(tabla_I);

disp(’Resultados para J_n (Método Estable):’);

disp(tabla_J);

ˆ El cálculo de In mediante recurrencia es estable solo para n ≤ 25. Para n > 25, se recomienda usar
integración numérica directa.
ˆ El cálculo de Jn por recurrencia no es fiable para ningún rango de n. Es preferible calcular Jn
directamente mediante la definición integral.

ˆ La solución propuesta en el código garantiza resultados estables y precisos para ambos casos.

Otro Método Alternativo para el Cálculo de In y Jn

Método Recursivo Descendente
Relación de recurrencia para In
Dado que la relación de recurrencia ascendente para In es:
n
In+1 = − (1 + (n − 1)In−1 ) ,
π2
podemos reescribirla para calcular In−1 en términos de In+1 y In :

−π 2 In+1 − n
In−1 = .
(n − 1)(n − 2)

Usaremos esta fórmula para calcular In de forma descendente, partiendo de los valores conocidos I40 e
I39 obtenidos mediante integración numérica.

27
Relación de recurrencia para Jn
De forma similar, la relación de recurrencia ascendente para Jn :

n(n + 1)
Jn+1 = Jn−1 ,
(n + 1)2 + 1

puede reescribirse como:

(n + 1)2 + 1
Jn−1 = Jn+1 .
n(n + 1)
Usaremos esta fórmula recursiva descendente, partiendo de los valores conocidos J40 y J39 obtenidos
mediante integración numérica.

Código MATLAB
El siguiente código implementa el cálculo recursivo descendente para In y Jn y compara los resultados
con los obtenidos mediante la función integral de MATLAB:
Las tablas generadas con el código incluirán los siguientes datos:
ˆ Valores de n de 40 a 0.
ˆ Valores calculados con el método iterativo para In y Jn .
ˆ Valores calculados mediante integración directa.
ˆ La diferencia absoluta entre los dos métodos.

% Código MATLAB para calcular I_n y J_n usando el método recursivo descendente
clear; clc;

% Constante piˆ2
pi2 = piˆ2;

% Valores iniciales para I_40 e I_39 usando la función integral

I40 = integral(@(x) x.ˆ40 .* cos(pi * x), 0, 1);
I39 = integral(@(x) x.ˆ39 .* cos(pi * x), 0, 1);

% Valores iniciales para J_40 y J_39 usando la función integral

J40 = integral(@(x) exp(-x) .* sin(x).ˆ40, 0, inf);
J39 = integral(@(x) exp(-x) .* sin(x).ˆ39, 0, inf);

% Inicialización de los vectores para almacenar resultados

I = zeros(1, 41); % I_n desde n=40 a n=0
J = zeros(1, 41); % J_n desde n=40 a n=0

% Asignar los valores iniciales conocidos

I(41) = I40; % I_40 en la posición 41 (ı́ndice MATLAB)
I(40) = I39; % I_39 en la posición 40
J(41) = J40; % J_40 en la posición 41
J(40) = J39; % J_39 en la posición 40

% Cálculo recursivo descendente para I_n

for n = 39:-1:1
I(n) = (-pi2 * I(n+2) - n) / (n * (n-1));
end

% Cálculo recursivo descendente para J_n

for n = 39:-1:1
J(n) = ((n+1)ˆ2 + 1) / (n * (n+1)) * J(n+2);
end

28
% Cálculo de I_n y J_n mediante la función integral
I_integral = zeros(1, 41);
J_integral = zeros(1, 41);

for n = 0:40
I_integral(n+1) = integral(@(x) x.ˆn .* cos(pi * x), 0, 1);
J_integral(n+1) = integral(@(x) exp(-x) .* sin(x).ˆn, 0, inf);
end

% Calcular diferencias
I_diff = abs(I - I_integral);
J_diff = abs(J - J_integral);

% Crear tablas de resultados

n_values = 40:-1:0;

I_table = table(n_values’, I’, I_integral’, I_diff’, ...

’VariableNames’, {’n’, ’I_n_Recursivo’, ’I_n_Integral’, ’Diferencia’});
J_table = table(n_values’, J’, J_integral’, J_diff’, ...
’VariableNames’, {’n’, ’J_n_Recursivo’, ’J_n_Integral’, ’Diferencia’});

% Mostrar resultados
disp(’Resultados para I_n:’);
disp(I_table);

disp(’Resultados para J_n:’);

disp(J_table);

Interpretación de resultados
El método recursivo descendente no es numéricamente estable para calcular . Sin embargo, los valores
calculados mediante el método iterativo coinciden con los obtenidos por integración directa dentro de
una tolerancia razonable, indicando la estabilidad y precisión del método recursivo descendente para .
Cualquier divergencia significativa se debe a errores de redondeo acumulados.

29