0% encontró este documento útil (0 votos)
7 vistas134 páginas

Métodos Numéricos en Ingeniería Civil

El documento aborda los métodos numéricos, su historia, y su importancia en la resolución de problemas matemáticos complejos mediante aproximaciones. Se discuten los algoritmos, errores en computación, y diversas técnicas para la solución de ecuaciones, interpolación, y cálculo integral y diferencial. Se enfatiza la necesidad de que los ingenieros comprendan y apliquen estos métodos en su práctica profesional.

Cargado por

marlonvalfren
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd
0% encontró este documento útil (0 votos)
7 vistas134 páginas

Métodos Numéricos en Ingeniería Civil

El documento aborda los métodos numéricos, su historia, y su importancia en la resolución de problemas matemáticos complejos mediante aproximaciones. Se discuten los algoritmos, errores en computación, y diversas técnicas para la solución de ecuaciones, interpolación, y cálculo integral y diferencial. Se enfatiza la necesidad de que los ingenieros comprendan y apliquen estos métodos en su práctica profesional.

Cargado por

marlonvalfren
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd

UNIVERSIDAD VERACRUZANA

FACULTAD DE INGENIERÍA CIVIL

EXPERIENCIA EDUCATIVA

METODOS NUMÉRICOS

MIC. VICTOR FERNÁNDEZ ROSALES

Última modificación: abril de 2021


Métodos numéricos

ÍNDICE

Página
INTRODUCCIÓN .............................................................................................................................. 1
Métodos numéricos ......................................................................................................................... 2
Historia ............................................................................................................................................ 3
Razones por las cuales se estudian los métodos numéricos ............................................................ 3
¿Dónde se utilizan? ......................................................................................................................... 4
Algoritmo ........................................................................................................................................ 4
Propiedades que deben cumplir los algoritmos numéricos ............................................................. 6
Convergencia............................................................................................................................... 6
Estabilidad ................................................................................................................................... 6
ERRORES ........................................................................................................................................... 7
Sistemas numéricos ......................................................................................................................... 7
Conversión de un número binario al sistema decimal ................................................................. 7
Manejo de números en la computadora........................................................................................... 8
Números enteros.......................................................................................................................... 8
Números reales (punto flotante) .................................................................................................. 9
Causas de errores graves en computación ..................................................................................... 10
Suma de números muy distintos en magnitud ........................................................................... 10
Resta de números casi iguales ................................................................................................... 10
Overflow y Underflow .............................................................................................................. 10
Tipos de errores ............................................................................................................................. 12
Error inherente........................................................................................................................... 12
Error de redondeo ...................................................................................................................... 12
Error por truncamiento .............................................................................................................. 12
Estimación del error por métodos iterativos.................................................................................. 14
Ejercicios ................................................................................................................................... 17
SOLUCIÓN DE ECUACIONES ALGEBRAICAS Y TRASCENDENTES ................................... 18
Método gráfico .............................................................................................................................. 19
Tipos de métodos numéricos para encontrar la raíz de una función ............................................. 21
Métodos cerrados o acotados ........................................................................................................ 21
Método de la bisección o método de bipartición del intervalo .................................................. 22

i
Métodos numéricos

Ejercicios ................................................................................................................................... 26
Método de la regla falsa o de la falsa posición.......................................................................... 28
Ejercicios ................................................................................................................................... 31
Métodos abiertos ........................................................................................................................... 32
Método del punto fijo ................................................................................................................ 32
Ejercicios ................................................................................................................................... 36
Método de Newton-Raphson ..................................................................................................... 37
Método de Newton-Raphson modificado para el cálculo de raíces múltiples .......................... 40
Ejercicios ................................................................................................................................... 42
Método de la Secante .................................................................................................................... 43
Ejercicios ................................................................................................................................... 46
Método de Müller .......................................................................................................................... 47
SOLUCIÓN DE ECUACIONES ALGEBRAICAS LINEALES SIMULTÁNEAS ........................ 50
Métodos directos ........................................................................................................................... 52
Método de eliminación de Gauss .............................................................................................. 52
Ejercicios ................................................................................................................................... 55
Método de Gauss-Jordan ............................................................................................................... 57
Ejercicios ................................................................................................................................... 61
Métodos Iterativos ......................................................................................................................... 62
Matriz diagonalmente dominante .............................................................................................. 62
Método de Jacobi o de los desplazamientos simultáneos .......................................................... 62
Ejercicios ................................................................................................................................... 67
Método de Gauss –Seidel .......................................................................................................... 69
Ejercicios ................................................................................................................................... 72
INTERPOLACIÓN ........................................................................................................................... 75
Aproximación polinomial simple e interpolación ......................................................................... 76
Ejercicios ................................................................................................................................... 80
Polinomios de Lagrange ................................................................................................................ 81
Ejercicios ................................................................................................................................... 85
Polinomio de Newton (Diferencias divididas) .............................................................................. 86
Ejercicios ................................................................................................................................... 90
INTEGRACIÓN NUMÉRICA ......................................................................................................... 91

ii
Métodos numéricos

Fórmulas de integración de Newton-Cotes ................................................................................... 92


Regla del Trapecio ........................................................................................................................ 93
Ejercicios ................................................................................................................................... 95
Regla del Trapecio compuesto ...................................................................................................... 96
Ejercicios ................................................................................................................................... 99
Regla de Simpson ........................................................................................................................ 100
Regla de Simpson 1/3.................................................................................................................. 100
Ejercicios ................................................................................................................................. 102
Regla de Simpson 1/3 Compuesta ............................................................................................... 103
Ejercicios ................................................................................................................................. 105
Regla de Simpson 3/8.................................................................................................................. 107
Ejercicios ................................................................................................................................. 109
Regla de Simpson 3/8 compuesta................................................................................................ 110
Ejercicios ................................................................................................................................. 113
Método de Romberg .................................................................................................................... 115
Ejercicios ................................................................................................................................. 118
SOLUCIÓN NUMÉRICA DE ECUACIONES DIFERENCIALES ORDINARIAS .................... 119
Método de Euler-Cauchy ............................................................................................................ 120
Ejercicios ................................................................................................................................. 123
Método de Runge-Kutta (Cuarto orden) ..................................................................................... 124
Ejercicios ................................................................................................................................. 127
REFERENCIAS .............................................................................................................................. 129

INDICE DE TABLAS

INDICE DE ECUACIONES

INDICE DE FIGURAS

Figura 1. Regla de cálculo (Elaboración propia) ................................................................................. 3


Figura 2. Representación de grafica de un algoritmo .......................................................................... 5
Figura 3. Graficas que muestran algunas formas de convergencia y divergencia ............................... 7

iii
Métodos numéricos

Figura 4. Representación gráfica de ecuación lineal con 2 y 3 incógnitas ........................................ 50


Figura 5. Representación gráfica de los tipos de sistemas de ecuaciones lineales (Kolman y Hill,
2006). ................................................................................................................................................ 51
Figura 6. Interpolación de la acetona con 2 nodos (Nieves y Domínguez, 2006) ............................. 77
Figura 7. Interpolación de la acetona con 3 nodos (Nieves y Domínguez, 2006) ............................. 78
Figura 8. Esquema del método de Euler.......................................................................................... 120

iv
INTRODUCCIÓN
Al momento de aplicar las matemáticas a situaciones del mundo real nos encontramos a
menudo con problemas que no pueden ser resueltos analíticamente o de manera exacta y
cuya solución debe ser abordada con ayuda de algún procedimiento numérico. A
continuación, consideramos algunos problemas típicos, ya formulados matemáticamente,
para los cuales estudiaremos técnicas numéricas de solución.

Problema 1. Encontrar el área de la región comprendida entre las gráficas de: 𝑦 = 2𝑠𝑒𝑛𝑥,
𝑦 = 𝑒 − 𝑥 con x ∈ [0, π].

Problema 2. Encontrar las raíces de la ecuación polinómica:


𝑥 5 + 11𝑥 4 − 21𝑥 3 − 10𝑥 2 − 21𝑥 − 5 = 0

Problema 3. Resolver los siguientes sistemas de ecuaciones:


a) El sistema lineal AX=b con:
2 −1 0 0 0 3
−1 2 −1 0 0 −2
𝐴 = 0 −1 2 −1 0 𝐵= 2
0 0 −1 2 −1 −2
[0 0 0 −1 2 ] [1]
b) El sistema no-lineal
𝑥 2 + 𝑥𝑦 3 = 9
{ 2
3𝑥 𝑦 − 𝑦 3 = 4

Problema 4. Dada la siguiente tabla de datos correspondiente a una cierta función


𝑦 = 𝑓(𝑥)
Nodo 1 2 3 4 5 6
xk -2 -1 0 1 2 3
f(xk) -5 1 1 1 7 25

Problema 5. Hallar el valor de cada una de las siguientes integrales:


𝜋
1 1 2 3
𝑠𝑒𝑛 𝑥 𝑥2
𝑠𝑒𝑛2 𝑥 1
∫ 𝑑𝑥 ∫ 𝑒 𝑑𝑥 ∫ √1 − 𝑑𝑥 (𝑒𝑙𝑖𝑝𝑡𝑖𝑐𝑎) ∫ 𝑑𝑥
𝑥 4 𝑙𝑛𝑥
0 0 0 2

En relación con los problemas anteriores, tenemos que:

1
Métodos numéricos

En el problema 1, es necesario determinar los puntos de intersección de las gráficas de 𝑦 =


2𝑠𝑒𝑛𝑥 y 𝑦 = 𝑒 − 𝑥, para lo se debe resolver la ecuación 2𝑠𝑒𝑛𝑥 = 𝑒 − 𝑥 y no se dispone
de un método algebraico para hacerlo.

En el problema 2, se trata de hallar los ceros de un polinomio de 5 grado y, como se sabe,


sólo se conocen métodos algebraicos para encontrar raíces de ecuaciones polinómicas de
grado menor o igual que 4.

En el problema 3, se tienen dos sistemas de ecuaciones: El de la parte a) es lineal y se


conocen métodos de solución (por ejemplo, el método de eliminación Gaussiana), sin
embargo, para sistemas de tamaño mayor, no sólo es conveniente sino necesario
implementar tales métodos a través del computador (método numérico). En la parte b) se
tiene un sistema no-lineal y no se conocen métodos algebraicos generales para resolverlo.

El problema 4, se puede resolver analíticamente (por interpolación), sin embargo, para


determinar los coeficientes de dichos polinomios existen técnicas que permiten encontrarlos
rápidamente y que pueden implementarse en el computador.

El problema 5, corresponde a integrales definidas cuyo integrando tiene antiderivada


(operación inversa a la derivada) que no es elemental (no existe ninguna fórmula en
términos de funciones elementales (es decir: polinomios, funciones trigonométricas,
exponenciales, logarítmicas y productos y composiciones de estas funciones)).

Métodos numéricos
Es un procedimiento mediante el cual se obtiene, casi siempre de manera aproximada, la
solución de ciertos problemas realizando cálculos puramente aritméticos y lógicos
(operaciones aritméticas elementales, cálculo de funciones, consulta de una tabla de
valores, cálculo proposicional, entre otros). Un tal procedimiento consiste en una lista
finita de instrucciones precisas que especifican una secuencia de operaciones algebraicas y
lógicas (algoritmo), que producen o bien una aproximación de la solución del problema
(solución numérica) o bien un mensaje. La eficiencia en el cálculo de dicha aproximación
depende, en parte, de la facilidad de implementación del algoritmo y de las características

2
Métodos numéricos

especiales y limitaciones de los instrumentos de cálculo (los computadores). En general, al


emplear estos instrumentos de cálculo se introducen errores llamados de redondeo.

Historia
Antes del uso o aparición de la PC, había 3 métodos diferentes que se aplican a la solución
de problemas:
1. Usando métodos exactos o analíticos (éstos tienen un valor práctico limitado ya que
son aplicables a una clase limitada de problemas).
2. Para analizar el comportamiento de los sistemas se usaban soluciones gráficas
(resultados no muy precisos, tediosos y difíciles de implementar sin ayuda de una
PC).
3. Para implementar los métodos numéricos se utilizaban calculadoras manuales y
reglas de cálculo como se muestra en la Figura 1, (son tediosos, lentos y no existen
resultados consistentes). Antes de la aparición y uso del PC se gastaba mucha
energía en la técnica misma de solución, en lugar de aplicarla sobre la definición del
problema y su interpretación.

Figura 1. Regla de cálculo (Elaboración propia)

Objetivo de su utilización
El objetivo principal del análisis numérico es encontrar soluciones “aproximadas” a
problemas complejos utilizando sólo las operaciones más simples de la aritmética. Se
requiere de una secuencia de operaciones algebraicas y lógicas que producen la
aproximación al problema matemático.

Razones por las cuales se estudian los métodos numéricos


• Es importante que el futuro ingeniero tenga los conocimientos básicos de los métodos
numéricos más comunes, ya que, en el transcurso de su carrera, tendrá la necesidad de

3
Métodos numéricos

usar software comercial o implementar su propio software, que resuelvan los


algoritmos de problemas reales y que estén basados sobre algún método numérico.
• Con los métodos numéricos el ingeniero usará la computadora como herramienta, el
cual es uno de los propósitos, porque el profesionista debe de olvidarse de los cálculos,
y enfocarse en el diseño y planteamiento de la solución de los problemas.
• Proporciona una mayor comprensión de las matemáticas, ya que reducen las
matemáticas superiores a operaciones básicas simples.

¿Dónde se utilizan?
Los métodos numéricos pueden ser aplicados para resolver procedimientos matemáticos en:
• Cálculo de derivadas
• Integrales
• Ecuaciones diferenciales
• Operaciones con matrices
• Interpolaciones
• Ajuste de curvas
• Polinomios, entre otros.

Algoritmo
Es una secuencia lógica de pasos necesarios para ejecutar una tarea específica tal como la
solución de un problema.

Características de los algoritmos


• Preciso: Definir de manera rigurosa, sin dar lugar a ambigüedades.
• Definido: Si se sigue un algoritmo dos veces, se obtendrá el mismo resultado.
• Finito: Debe terminar en algún momento.
• Puede tener cero o más elementos de entrada.
• Debe producir un resultado: Los datos de salida serán los resultados de efectuar las
instrucciones.

Etapas para la solución de un problema


1. Análisis del problema, definición y delimitación (macro algoritmo). Considerar los
datos de entrada, el proceso que debe realizar la computadora y los datos de salida.

4
Métodos numéricos

2. Diseño y desarrollo del algoritmo (se utiliza seudocódigo, escritura natural del
algoritmo, diagramas de flujo, entre otros).
3. Prueba de escritorio. Seguimiento manual de los pasos descritos en el algoritmo. Se
hace con valores bajos y tiene como fin detectar errores.
4. Codificación. Selección de un lenguaje de programación y digitalización del
seudocódigo haciendo uso de la sintaxis y estructura gramatical del lenguaje
seleccionado.
5. Compilación o interpretación del programa. El software elegido convierte las
instrucciones escritas en el lenguaje a las comprendidas por la computadora.
6. Ejecución. El programa es ejecutado por la computadora para llegar a los resultados
esperados.
7. Depuración. (debug). Operación de detectar, localizar y eliminar errores de mal
funcionamiento del programa.
8. Evaluación de resultados. Obtenidos los resultados se los evalúa para verificar si
son correctos. (Un programa puede arrojar resultados incorrectos aun cuando su
ejecución no muestre errores)

Un algoritmo se puede representar mediante un diagrama de flujo la cual nos da una


representación gráfica en la cual se emplean bloques y flechas.

Diseño del Programa de


Problema algoritmo computadora

Figura 2. Representación de grafica de un algoritmo


Ejemplo:
Algoritmo para la solución de la suma de 2 números cualquiera.
1. Inicio
2. Solicitar el valor de a
3. Solicitar el valor de b
4. Sumar a con b y asignar a c la respuesta
5. Imprimir el valor de c
6. fin

5
Métodos numéricos

Propiedades que deben cumplir los algoritmos numéricos

Convergencia
Se entiende por convergencia la garantía de que, al realizar un número de repeticiones
(iteraciones), las aproximaciones obtenidas terminan por acercarse cada vez más al valor
verdadero como se muestra en la Figura 3.

Estabilidad
Cuando en un algoritmo o método numérico el crecimiento de los errores que traen los
datos es lineal, entonces se dice que el algoritmo o método numérico es estable y los
resultados que nos arroje serán válidos, por el contrario, si el crecimiento del error es
exponencial entonces el algoritmo es inestable y no puede tomarse como válidos los
resultados obtenidos.

6
Métodos numéricos

Figura 3. Graficas que muestran algunas formas de convergencia y divergencia

ERRORES
Los métodos numéricos ofrecen soluciones aproximadas muy cercanas a las soluciones
exactas; la discrepancia entre una solución verdadera y una aproximada constituye un
error, por lo que es importante saber qué se entiende por aproximar y aprender a cuantificar
los errores para minimizarlos.

Sistemas numéricos
Conversión de un número binario al sistema decimal
Teniendo en cuenta el valor de cada dígito en su posición, que es el de una potencia de 2,
cuyo exponente es 0 en el bit situado más a la derecha, y se incrementa en una unidad
según vamos avanzando posiciones hacia la izquierda.
Ejemplos:
10110012= 1010102=
1x26 + 0x25+ 1x24+ 1x23+ 0x22+ 0x21+ 1x20= 1x25+ 0x24+ 1x23+ 0x22+ 1x21+ 0x20=
64 + 0 + 16 + 8 + 0 + 0 + 1 =8910 32 + 0 + 8 + 0 + 2 + 0 =4210

Conversión de números enteros del sistema decimal al sistema binario


De una palabra de 16 bits el primero debe ser el signo donde (0) es positivo y (1) es
negativo y los 15 bits restantes son para el número, esto nos da un intervalo de:
0000000000000000 a 1111111111111111 que en decimal es de -32,767 a 32,767.

7
Métodos numéricos

Método: Se realizan divisiones sucesivas por 2 y se escriben los restos obtenidos en cada
división en orden inverso al que han sido obtenidos. (Nota: el residuo se multiplica por 2
para sacar el resto).
Ejemplos: 4710= (en una palabra, de 16 bits) =1011112
47/2 =23 resto 1
23/2 =11 resto 1
11/2 =5 resto 1
5/2 =2 resto 1
2/2 = 1 resto 0
1 / 2 = 0 resto 1
0 0 0 0 0 0 0 0 0 0 1 0 1 1 1 1

52510= (en una palabra, de 16 bits) =10000011012


525/2 =262 resto 1
262/2 =131 resto 0
131/2 =65 resto 1
65/2 =32 resto 1
32/2 = 16 resto 0
16/2 = 8 resto 0
8/2 = 4 resto 0
4/2 = 2 resto 0
2/2 = 1 resto 0
1/2 = 0 resto 1
0 0 0 0 0 0 1 0 0 0 0 0 1 1 0 1

Manejo de números en la computadora


Por razones prácticas, solo se puede manejarse una cantidad finita de bits para cada número
en una computadora, y esa cantidad o longitud varía de una maquina a otra.

Para una computadora dada, el número de bits generalmente se llama palabra. Las palabras
van de 8 hasta 64 bits, por ejemplo, una palabra de 32 bits puede dividirse en 4 bytes (8
bites cada una).

Números enteros
Cada palabra, cualquiera que sea su longitud, almacena un número, aunque en ciertas
circunstancias se usan varias palabras para contener un número. Por ejemplo, considérese
una palabra de 16 bits para almacenar números enteros. De los 16 bits, el primero
representa el signo del número; un cero es un signo más y un uno un signo menos. Los 15

8
Métodos numéricos

bits restantes pueden usarse para guardar números binarios desde 0000000000000000 hasta
1111111111111111.

Ejemplo:
-12510= -11111012

1 0 0 0 0 0 0 0 0 1 1 1 1 1 0 1

Números reales (punto flotante)


Para almacenar un número real se emplea en su representación binaria, llamada de punto
flotante la notación.

En los bits del 1 al 7 se almacena el exponente de la base 2 y en los 8 bits restantes la


fracción.

Según en el lenguaje de los logaritmos, la fracción es llamada mantisa y la exponente


característica.

Ejemplo:
-125.3210= -1111101.0101000111101012

Normalizado queda (en base a la especificación de IEEE 754):

-.1111101010100011110101 X 2 + 111 (7 en binario, que fue el corrimiento del punto a la izquierda)

bits truncados en el almacenamiento

Signo mantisa (número)

Característica positiva (exponente)

1 0 0 0 0 1 1 1 1 1 1 1 1 0 1 0
Característica Mantisa

9
Métodos numéricos

Causas de errores graves en computación

Suma de números muy distintos en magnitud

Ejemplo:
0.002 =0.2000 x 10-2
600 =0.6000 x 103

Números normalizados
0.000002 x 103
+ 0.600000 x 103
0.600002 x 103

Como solo se puede manejar 4 dígitos, los últimos 2 son eliminados y la suma es 0.6000 x
103, por lo que la suma nunca se realizó.

Resta de números casi iguales


Ejemplo:
.2145 x 100
- .2144 x 100
.0001 x 100
Que normalizado el resultado es: 0.1 x 10-3

Como solo existe un digito significativo se sugiere no confiar en su exactitud.

Overflow y Underflow
Con frecuencia una operación aritmética con 2 números válidos da como resultado un
número tan grande o pequeño que la computadora no puede manejarlo.
Ejemplo:
0.5000 x 108
x 0.2000x 109
0.1000 x 1017
El producto es muy grande y no puede almacenarse porque la característica requiere 3
dígitos y se produce un overflow.

El underflow puede aparecer en la multiplicación o división, y por lo general no es tan serio


como el overflow; las computadoras casi nunca envían mensajes de underflow.

10
Métodos numéricos

Ejemplo:
Cuando se suma 10,000 veces 0.0001 con él mismo, debe resultar 1; sin embargo, el
número 0.0001 en binario resulta en una sucesión infinita de ceros y unos que se trunca al
ser almacenada en una palabra de memoria, con lo que se perderá información y el
resultado de la suma ya no será 1.
function error=error()
% function para demostrar el error
format long
s=0;
for i=1:10000
s=s+0.0001;
end
disp(s)

s=1;
for i=1:10000
s=s+0.0001;
end
disp(s)

s=1000;
for i=1:10000
s=s+0.0001;
end
disp(s)

s=10000;
for i=1:10000
s=s+0.0001;
end
disp(s)
end

Resultados:

0.999999999999906
1.999999999999890
1.000999999999749e+003
1.000099999999293e+004

11
Métodos numéricos

Tipos de errores
Error inherente
En muchas ocasiones, los datos con que se inician los cálculos contienen un cierto error
debido a que se han obtenido mediante la medida experimental de una determinada
magnitud física. Así, por ejemplo, el diámetro de la sección de una varilla de acero
presentaría un error según se haya medido con una cinta métrica o con un pie de rey. A este
tipo de error se le denomina error inherente.

Error de redondeo
Como no es posible guardar un número binario de longitud infinita o un número de más
dígitos de los que posee la mantisa de la computadora que se está empleando, se almacena
sólo un número finito de estos dígitos; como consecuencia, se comete automáticamente un
pequeño error, conocido como error de redondeo, que al repetirse muchas veces puede
llegar a ser considerable.

Error por truncamiento


Los errores por truncamiento ocurren cuando un número, cuya parte fraccionaria está
constituida por un número infinito de dígitos, requiere ser representado numéricamente en
forma aproximada, utilizando un número de cifras significativas.
Por ejemplo 3.1416 es una buena aproximación del número π, pero el valor exacto no
puede ser expresado numéricamente por completo, pues consta de un número infinito de
dígitos: 3.1415926535…; lo mismo ocurre con el 2.7183 para el número e, el 1.4142 para
√2, y el 0.333333 para 1/3.
Sin embargo, todos los números, ya sean enteros, racionales o irracionales, pueden ser
representados a través de formulaciones matemáticas exactas, utilizando series infinitas;
obviamente, las representaciones numéricas acotadas a un determinado número de cifras
significativas son aproximaciones numéricas que llevan implícitos errores por
truncamiento.

Por ejemplo, los números 1, 1/3 y e pueden expresarse matemáticamente, de manera exacta,
a través de las siguientes series infinitas:

12
Métodos numéricos

1 1 1 1 1 1
1= + + + + + +⋯
2 4 8 16 32 64
1 3 3 3 3 3
= + + + + +⋯
3 10 100 1000 10000 1000000
1 1 1 1 1 1
𝑒= + + + + + +⋯
0! 1! 2! 3! 4! 5!

Para ambos tipos de errores la relación entre el resultado exacto, o verdadero y el


aproximado está dado por:

valor verdadero = valor aproximado + error

El error numérico o error verdadero se define como:


𝐸𝑡 = 𝑣𝑎𝑙𝑜𝑟 𝑣𝑒𝑟𝑑𝑎𝑑𝑒𝑟𝑜 − 𝑣𝑎𝑙𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 [𝐸𝑟𝑟𝑜𝑟 𝑛𝑢𝑚é𝑟𝑖𝑐𝑜]

El error relativo porcentual se define como:


𝑒𝑟𝑟𝑜𝑟 𝑣𝑒𝑟𝑑𝑎𝑑𝑒𝑟𝑜
𝐸𝑟 = 100% [𝐸𝑟𝑟𝑜𝑟 𝑟𝑒𝑙𝑎𝑡𝑖𝑣𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑣𝑎𝑙𝑜𝑟 𝑣𝑒𝑟𝑑𝑎𝑑𝑒𝑟𝑜

Ejemplo:
Suponga que se tiene que medir la longitud de un puente y la de un remache, y se obtiene
9,999 y 9 cm respectivamente. Si los valores verdaderos son 10,000 y 10 cm, calcule:

a) el error numérico
b) el error relativo porcentual en cada caso:

Solución:
a) El error numérico en la medición del puente es: Et= 10,000 - 9,999 = 1 cm
y en la del remache es de: Et = 10 - 9 =1 cm.

b) El error relativo porcentual para el puente es:


1
𝐸𝑟 = 100% = 0.01%
10,000

y para el remache es de:


1
𝐸𝑟 = 100% = 10%
10
Por tanto, aunque ambas medidas tienen un error numérico de 1 cm, el error relativo
porcentual del remache es mucho mayor. Se concluye entonces que se ha hecho un buen
trabajo en la medición del puente; mientras que la estimación para el remache dejó mucho
que desear.

13
Métodos numéricos

Estimación del error por métodos iterativos


Uno de los retos que enfrentan los métodos numéricos es el de determinar estimaciones del
error en ausencia del conocimiento de los valores verdaderos. Por ejemplo, ciertos métodos
numéricos usan un método iterativo para calcular los resultados. En tales métodos se hace
una aproximación considerando la aproximación anterior. Este proceso se efectúa varias
veces, o de forma iterativa, para calcular en forma sucesiva, esperando cada vez mejores
aproximaciones. En tales casos, el error a menudo se calcula como la diferencia entre la
aproximación previa y la actual. Por tanto, el error aproximado porcentual está dado por:

𝒂𝒑𝒓𝒐𝒙𝒊𝒎𝒂𝒄𝒊ó𝒏 𝒂𝒄𝒕𝒖𝒂𝒍 − 𝒂𝒑𝒓𝒐𝒙𝒊𝒎𝒂𝒄𝒊ó𝒏 𝒂𝒏𝒕𝒆𝒓𝒊𝒐𝒓


𝑬𝒂 = | | 𝟏𝟎𝟎 % [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝒂𝒑𝒓𝒐𝒙𝒊𝒎𝒂𝒄𝒊ó𝒏 𝒂𝒄𝒕𝒖𝒂𝒍

Un caso muy interesante es una investigación que realiza Scarborough 1966, en que
determinó el número de cifras significativas que contiene el error como: Si reemplazamos
el Error esperado (Es) en la ecuación obtendremos el número de cifras significativas en
que es confiable el valor aproximado obtenido.

𝑬𝒔 = (𝟎. 𝟓 𝒙 𝟏𝟎𝟐−𝒏 )% 𝑛 = 𝑛ú𝑚𝑒𝑟𝑜 𝑑𝑒 𝑑𝑒𝑐𝑖𝑚𝑎𝑙𝑒𝑠 [𝐸𝑟𝑟𝑜𝑟 𝑒𝑠𝑝𝑒𝑟𝑎𝑑𝑜]

Planteamiento del problema: En matemáticas con frecuencia las funciones se representan


mediante series infinitas. Por ejemplo, el cálculo exacto de la función exponencial de un
número mediante la Serie de Taylor requiere infinitos sumandos.

𝑥
𝑥2 𝑥3 𝑥𝑛
𝑒 = 1 + 𝑥 + + +. . +
2! 3! 𝑛!

Así cuantos más términos se le agreguen a la serie, la aproximación será cada vez más una
mejor estimación del valor verdadero de ex.

Ejemplo:
Si calculamos la serie anterior con el valor de e con x=0.5 para la función con un error
menor al 0.05.

Calculando el error numérico y el error aproximado porcentual (Ea) empezando con el 1er
termino y posteriormente agregando más términos (1º, 2º, 3º, entre otros) hasta que el valor

14
Métodos numéricos

absoluto del valor aproximado sea menor al criterio prestablecido, que es que contemple 5
cifras significativas.

Solución:
𝑒 0.5 = 1.648721

En primer lugar, la ecuación del error esperado (Es) se emplea para determinar el criterio
de error que asegura que un resultado sea correcto en al menos tres cifras significativas:
𝐸𝑠 = (0.5 𝑥 10 2−5 )% = 0.5 𝑥 10−3 = 0.0005% = 0.05
para x=0.5

1er termino: 𝑒 0.5 = 1 =1


(1−0)
𝐸𝑎 = | | ∗ 100% = 1% = 100
1

2o termino: 𝑒 0.5 = 1 + 𝑥 =
𝑒 0.5 = 1 + 0.5 = 𝟏. 𝟓
(1.5−1)
𝐸𝑎 = | 1.5 | ∗ 100% = 0.3333% = 33.333

𝑥2
3er termino: 𝑒 0.5 = 1 + 𝑥 + =
2!
(0.5)2
𝑒 0.5 = 1 + 0.5 + = 𝟏. 𝟔𝟐𝟓
2!
(1.625−1.5)
𝐸𝑎 = | | ∗ 100% = 0.0769% = 7.692
1.625

𝑥2 𝑥3
4o termino: 𝑒 0.5 = 1 + 𝑥 + + =
2! 3!
(0.5)2 (0.5)3
𝑒 0.5 = 1 + 0.5 + + = 𝟏. 𝟔𝟒𝟓𝟖𝟑𝟑𝟑𝟑𝟑
2! 3!
(1.645833333−1.625)
𝐸𝑎 = | | ∗ 100% = 0.0127% = 1.266
1.645833333

𝑥2 𝑥3 𝒙𝟒
5o termino: 𝑒 0.5 = 1 + 𝑥 + + + =
2! 3! 4!
(0.5)2 (0.5)3 (0.5)4
𝑒 0.5 = 1 + 0.5 + + + = 𝟏. 𝟔𝟒𝟖𝟒𝟑𝟕𝟓𝟎𝟎
2! 3! 4!
(1.648437500−1.645833333)
𝐸𝑎 = | | ∗ 100% = 0.00158% = 0.158
1.648437500

𝑥2 𝑥3 𝒙𝟒 𝒙𝟓
6o termino: 𝑒 0.5 = 1 + 𝑥 + + + + =
2! 3! 4! 5!
(0.5)2 (0.5)3 (0.5)4 (0.5)5
𝑒 0.5 = 1 + 0.5 + + + + = 𝟏. 𝟔𝟒𝟖𝟔𝟗𝟕𝟗𝟏𝟕
2! 3! 4! 5!
(1.648697917−1.648437500)
𝐸𝑎 = | 1.648697917
| ∗ 100% = 0.000158% = 0.016

15
Métodos numéricos

Términos Resultado Ea (%) Es(%)


1 1 100 0.05
2 1.5 33.3
3 1.625 7.69
4 1.645833333 1.27
5 1.648437500 0.158
6 1.648697917 0.0158 parar

Así después de usar seis términos, el error aproximado porcentual (Ea) =0.0158 es menor
que el error esperado (Es)=0.05 y el cálculo termina.

16
Métodos numéricos

Ejercicios
1. Convierta los siguientes números decimales a los sistemas con base 2.
a) 536 b) 923 c) 1536

2. Convierta los siguientes números dados en binario a sistema decimal.


a) 1000 b)10101 c)111111

3. Calcule el valor de e a la -8.3 para la función:


−𝑥
𝑥2 𝑥3
𝑒 = 1 − 𝑥 + − +..
2! 3!
1 1
𝑒 −𝑥 = 𝑥
= 𝑥 2 𝑥3
𝑒 1+𝑥+ + +⋯
2! 3!

Compare con el valor verdadero de 4,023.87239 y comente los resultados. Use 25


términos para evaluar cada serie.

17
Métodos numéricos

SOLUCIÓN DE ECUACIONES ALGEBRAICAS Y TRASCENDENTES

Introducción
La determinación de las raíces de una ecuación es uno de los problemas más antiguos en
matemáticas y se han realizado un gran número de esfuerzos en este sentido. Su
importancia radica en que si podemos determinar las raíces de una ecuación también
podemos determinar máximos y mínimos, valores propios de matrices, resolver sistemas de
ecuaciones lineales y diferenciales, entre otros.

¿Qué es una raíz?


La raíz de una ecuación es aquel valor de la variable independiente que hace que el
resultado de la ecuación sea cero o por lo menos se acerque a cero con un cierto grado de
aproximación deseado (error máximo permitido).

El método numérico se utiliza en una función algebraica no lineal cuando no se puede


despejar la variable que interesa estudiar.

Antes de la llegada de las computadoras digitales se disponía de una serie de métodos para
encontrar las raíces de ecuaciones algebraicas y trascendentes. En algunos casos las raíces
−𝑏±√𝑏 2 −4𝑎𝑐
se obtenían con métodos directos, como se hace con la ecuación 𝑥 = para
2𝑎

resolver 𝑓(𝑥) = 𝑎𝑥 2 + 𝑏𝑥 + 𝑐 = 0.

Sin embargo, existen ecuaciones que no se resuelven directamente y aparecen muchas más
en las que no es posible encontrar solución. Por ejemplo, incluso una función tan simple
como 𝑓(𝑥) = 𝑒 −𝑥 – 𝑥 no se puede resolver en forma analítica. En tales casos, la única
alternativa es una técnica con solución aproximada.

Los polinomios son un caso simple


f(x)= a0 + a1x + a2x2 + … + anxn donde las a son constantes

Ejemplos:
𝑓(𝑥) = 1 − 2.5𝑥 + 7𝑥 2
𝑓(𝑥) = 3𝑥 2 − 𝑥 3 + 7𝑥 5

18
Métodos numéricos

Una función transcendental es una función que no es algebraica, incluye funciones


trigonométricas, exponenciales, logarítmicas entre otras, por ejemplo:

2
𝑓(𝑥) = 𝑙𝑛 𝑥 − 1, 𝑓(𝑥) = 𝑒 −0.2𝑥 𝑦 𝑓(𝑥) = 𝑆𝑒𝑛 (3𝑥 − 0.5)

1. La determinación de raíces reales de ecuaciones algebraicas y transcendentales:


Están diseñadas para determinar el valor de una raíz simple de acuerdo con un
conocimiento previo de su posición aproximada.
2. La determinación de todas las raíces reales y complejas de un polinomio: Están
diseñados específicamente para polinomios, determinan sistemáticamente todas las
raíces en lugar de una, dada una aproximación según una posición.

Un método para obtener una solución aproximada consiste en graficar la función para
determinar dónde cruza el eje de las x. Este punto, que representa el valor de x para el cual
𝑓(𝑥) = 0 es la raíz.

Método gráfico
Un método simple para obtener una aproximación a la raíz de la ecuación y observar en
donde cruza el eje X. Este punto que representa el valor de x para el cual 𝑓(𝑥) = 0
proporciona la aproximación inicial de la raíz.

Aunque los métodos gráficos son útiles en la obtención de estimaciones de las raíces, tienen
el inconveniente de que son poco precisos.
Otro método es el de prueba y error. Esta técnica consiste en elegir un valor de x y evaluar
si f(x) es cero. Si no es así se hace otra elección y se evalúa nuevamente f(x) para
determinar si el valor ofrece una mejor aproximación de la raíz. El proceso se repite hasta
que se obtenga un valor que proporcione una f(x) cercana a cero; por lo tanto, se crearon
métodos más exactos y fáciles de adoptarlos a las computadoras, reduciendo así el tiempo
en encontrar la solución y en la exactitud de estos.

Ejemplo:
Utilizar el método Gráfico para encontrar la raíz de la ecuación: 𝑓(𝑥) = 𝑒 −𝑥 − 𝑥

19
Métodos numéricos

x f(x) cálculos
-0.4 1.89 e0.4+ 0.4=1.89
-0.2 1.42 e0.2 + 0.2=1.42
0 1.00 e0 + 0=1.00
0.2 0.62 e-0.2 – 0.2=0.62
0.4 0.27 e-0.4 – 0.4=0.27
0.6 -0.05 e-0.6 – 0.6=-0.05

f(x)
2

1.5

0.5

0
-0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8
-0.5

La raíz se encuentra en: x= 0.56


Las técnicas graficas tienen un valor práctico limitado, ya que no son precisas, estas solo se
usarán para obtener las aproximaciones, las cuales se pueden ocupar como valores iniciales
en los métodos numéricos.
Las interpretaciones geométricas además de aproximar a la raíz son herramientas
importantes en el aislamiento de las propiedades de las funciones previendo las fallas de los
métodos numéricos, en general si f(x1) y f(xu) tienen signos opuestos existe un número
impar de raíces dentro del intervalo definido por los mismos, si f(xi) y f(xy) tienen el mismo
signo, no hay raíces o existe un número par de ellas entre los valores dados.
Por medio de las computadoras se puede acelerar las soluciones aproximadas y obtener sus
características, disminuyendo el tiempo y la exactitud que se si hacen de forma manual.

Ejercicio:
Graficar la función: 𝑓(𝑥) = 𝑠𝑒𝑛10𝑥 + 𝑐𝑜𝑠3𝑥 desde x = [-1:0.5:1] en una hoja de cálculo y
software matemático. ¿Cuántas raíces se encontraron?
Y si ahora se modifica el rango a x = [-1:0.2:1] ¿Cuantas raíces se encontraron?

20
Métodos numéricos

Tipos de métodos numéricos para encontrar la raíz de una función

𝐵𝑖𝑠𝑒𝑐𝑐𝑖ó𝑛
Cerradas o acotados: {
𝐹𝑎𝑙𝑠𝑎 𝑝𝑜𝑠𝑖𝑐𝑖ó𝑛
(Requieren de dos valores de
x que encierren a la raíz)

𝑃𝑢𝑛𝑡𝑜𝑓𝑖𝑗𝑜
Abiertos: {𝑁𝑒𝑤𝑡𝑜𝑛 − 𝑅𝑎𝑝ℎ𝑠𝑜𝑛
𝑆𝑒𝑐𝑎𝑛𝑡𝑒
𝑀𝑢𝑙𝑙𝑒𝑟
(Requieren de un valor de x o dos, pero
que no necesariamente encierren a la raíz)

La importancia de utilizar el método grafico en los métodos numéricos para encontrar


raíces de la función, es, ayudar a encontrar tanto las raíces reales como las imaginarias al
presentar el comportamiento de la función en un eje cartesiano.

Métodos cerrados o acotados


Son métodos que aprovechan el hecho de que una función en forma típica cambia de signo
en la cercanía de una raíz, ya que necesitan de 2 valores iniciales para la raíz. Como su
nombre indica, esto valores deben “encerrar” o estar sobre cualquier lado de la raíz.
Emplean diferentes estrategias para reducir sistemáticamente el tamaño del intervalo y así
converger a la respuesta correcta.

21
Métodos numéricos

Método de la bisección o método de bipartición del intervalo


Éste se clasifica como un método de acotamiento. Es aplicable a ecuaciones de la forma
𝑓(𝑥) = 0, cuando es posible encontrar dos valores limitantes xa y xb tales que la función
f(x) cambia de signo una vez para valores x en el intervalo (xa ≤ x ≤ xb). Por consiguiente,
los valores limitantes acotan la raíz.

En cada iteración el tamaño del intervalo se reduce a la mitad después de n iteraciones, el


intervalo original se habrá reducido 2𝑛 veces.

y f(x)
f(xb)

f(xr1)
xa=xa1 xr2 xr3 xr
x
xr1 Xb=xb
f(xr3)
1
f(xr2)
f(xa) xa1 xr1 xb1
xa2 xr2 xb2
xa3 xr3 xb3

𝒙𝒂 + 𝒙𝒃
𝒙𝒓 = ( ) [𝐵𝑖𝑠𝑒𝑐𝑐𝑖ó𝑛]
𝟐

El error aproximado porcentual está dado por:


𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑐𝑖ó𝑛 𝑎𝑐𝑡𝑢𝑎𝑙 − 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑐𝑖ó𝑛 𝑎𝑛𝑡𝑒𝑟𝑖𝑜𝑟
𝐸𝑎 = | | 100% [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑐𝑖ó𝑛 𝑎𝑐𝑡𝑢𝑎𝑙

22
Métodos numéricos

Algoritmo:
1. Elijase los valores iniciales xa y xb de tal manera que la función cambie de signo
sobre el intervalo, esto se puede verificar asegurándose de que: f(xa).f(xb)<0
2. La primera aproximación a la raíz se determina como:
(𝑥𝑎 + 𝑥𝑏 )
𝑥𝑟 =
2
3. Realice las siguientes evaluaciones y determine en que subintervalo se ubica la raíz.
• Si f(xa).f(xr)<0 entonces la raíz se encuentra dentro del primer sub intervalo, por
lo tanto resuélvase xb=xr y continuase con el paso 4.
• Si f(xa).f(xr) >0 entonces la raíz se encuentra en el segundo sub intervalo, por lo
tanto resuélvase xa=xr y continuase con el paso 4.
• Si f(xa).f(xr)=0 entonces la raíz es igual a xr y se terminan los cálculos.
4. Calcule una nueva aproximación a la raíz como en el paso 2.
5. Asegúrese que la nueva aproximación es tan exacta como se desea, si es así
entonces los cálculos terminan, de lo contrario regrese al paso 3.

Ejemplo:
Determinar por el método de la bisección la raíz de la siguiente función:
𝑓(𝑥) = 𝑒 −𝑥 − 𝑥, tomando 5 decimales.

𝐸𝑠 = (0.5 𝑥 102−𝑛 )%
Es = (0.5 x 102-5)% = (0.5 x 10-3)% = 0.0005%= 0.05

Método de sustitución
x f(x) resultado f(x)
-1 3.7183 exp(1) +1=3.7183
4
-0.5 2.1487 exp(0.5) +0.5=2.1487
0 1 exp(0) -0=1
0.5 0.1065 exp(-0.5) -0.5=0.1065 2
1 -0.6321 exp(-1) -1=-0.6321
1.5 -1.2769 exp(-1.5) -1.5=-1.2769 0
-1.5 -1 -0.5 0 0.5 1 1.5 2
-2

23
Métodos numéricos

f(xa).f(xr)<0 f(xa).f(xr)>0
iter xa xb xr f(xa) f(xr) f(xa).f(xr)=0 Ea(%) Es(%)
xb=xr xa=xr
1 0.5 1.0 0.75000 0.10653 -0.27763 si no no 100 0.05
2 0.50000 0.75000 0.62500 0.10653 -0.08974 si no no 20.00000
3 0.50000 0.62500 0.56250 0.10653 0.00728 no si no 11.11111
4 0.56250 0.62500 0.59375 0.00728 -0.04150 si no no 5.26316
5 0.56250 0.59375 0.57813 0.00728 -0.01718 si no no 2.70270
6 0.56250 0.57813 0.57031 0.00728 -0.00496 si no no 1.36986
7 0.56250 0.57031 0.56641 0.00728 0.00116 no si no 0.68966
8 0.56641 0.57031 0.56836 0.00116 -0.00191 si no no 0.34364
9 0.56641 0.56836 0.56738 0.00116 -0.00038 si no no 0.17212
10 0.56641 0.56738 0.56689 0.00116 0.00039 no si no 0.08613
11 0.56689 0.56738 0.56714 0.00039 0.00001 no si no 0.04305 parar

Iteración 1
xr= (0.5+1)/2 = 0.75
f(xa)=e-0.5 -0.5= 0.1065
f(xr)=e-0.75 – 0.75 = -0.2776
Ea =|(0.75-0)/0.75| % = 1% = 100

Iteración 2
xr= (0.5 + 0.75)/2 = 0.625
f(xa)=e-0.5 -0.5 = 0.1065
f(xr)=e-0.625 – 0.625 = -0.0.897
Ea =|(0.625 – 0.75)/0.625| % = 0.2000% = 20.00

Iteración 3
xr= (0.5 + 0.625)/2 = 0.5625
f(xa)=e-0.5 -0.5 = 0.1065
f(xr)=e-0.5625 – 0.5625 = 0.0073
Ea =|(0.5625 – 0.625)/0.5625| % =0.1111% = 11.11

Iteración 4
xr= (0.5625 + 0.625)/2 = 0.5938
f(xa)=e-0.5625 -0.5625 = 0.0073
f(xr)=e-0.5938 – 0.5938 = -0.0416
Ea =|(0.5938 - 0.5625)/0.5938| % =0.0527% = 5.27

Iteración 5
xr= (0.5625 + 0.5938)/2 = 0.5781
f(xa)=e-0.5625 -0.5625 = 0.0073
f(xr)=e-0.5781 – 0.5781 = -0.0171
Ea =|(0.5781 -0.5938)/0.5781| % = 0.0272% = 2.72

Nota: Se continúan las iteraciones hasta que Ea<= 0.05

24
Métodos numéricos

Algoritmo para el método de la bisección


function biseccion=biseccion()
% función para calcular la raíz por medio del método de bisección
clc;
f =input('Dame la función : ','s');
xa =input('Dame el intervalo inferior :');
xb =input('Dame el intervalo superior :');
err =input('Dame el porciento de error :');
f =inline(f);
if f(xa)*f(xb)<0
ezplot(f), grid on
% axis([xa xb -10 310]); % Asigna los límites de los ejes
% xlabel('pies'); % Asigna la etiqueta en el eje x
% ylabel('pies'); % Asigna la etiqueta en el eje y
erro=100;
xr =0;
i =0;
while erro>err
i=i+1;
ea=xr;
xr=(xa + xb)/2;
if f(xa)*f(xr)>0
xa=xr;
else
xb=xr;
end
erro=abs(((ea-xr)/xr)*100);
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',xr,i);
else
fprintf('\n\nNo existe la raíz en el intervalo\n')
end

Solución en Matlab:
Dame la función: exp(-x) -x
Dame el intervalo inferior: 0.5
Dame el intervalo superior: 1
Dame el porciento de error: 0.05

Resultado de la raíz: 0.567139 en 11 iteraciones


exp(-x)-x

200

150

100

50

-6 -4 -2 0 2 4 6
x

25
Métodos numéricos

Ejercicios
1. Determinar por el método de la bisección la raíz de las siguientes funciones:
1) 𝑓(𝑥) = 𝑥 3 + 2𝑥 2 + 10𝑥 − 20 tomando 5 decimales
2) 𝑓(𝑥) = −0.5𝑥 2 + 2.5𝑥 + 4.5 tomando 5 decimales
3) 𝑓(𝑥) = 𝑠𝑒𝑛(𝑥) + 0.8 𝑐𝑜𝑠(𝑥) en el intervalo [2,3] tomando 5 decimales
4) 𝑓(𝑥) = 𝑥 2 − 4𝑥 + 3.5 − 𝐼𝑛(𝑥) en [1,3] tomando 5 decimales
5) 𝑓(𝑥) = (𝑥 − 2.1)2 − 7𝑥𝑐𝑜𝑠(𝑥) en [1,2] tomando 5 decimales
6) 𝑓(𝑥) = 𝑠𝑒𝑛(𝑥) − 𝑐𝑜𝑠(1 + 𝑥 2 ) − 1 en [2π/3, π] tomando 5 decimales
7) 𝑓(𝑥) = 𝑥 2 𝑙𝑛𝑥 − 9𝑥 − 18 en el intervalo [6, 7] tomando 5 decimales
8) 𝑓(𝑥) = 𝑥 3 − 2𝑥 2 𝑠𝑒𝑛(𝑥) en el intervalo [2,3] tomando 5 decimales
9) 𝑓(𝑥) = 2𝑥 𝑡𝑎𝑛(𝑥) − 10 en el intervalo [1.3, 1.4] tomando 5 decimales
10) 𝑓(𝑥) = 𝑥 𝑙𝑜𝑔 𝑥 − 10 en el intervalo [5,6] tomando 5 decimales

1. Cuando se requiere encontrar la acidez de una solución de hidróxido de magnesio en ácido


clorhídrico, se obtiene la siguiente ecuación:

𝑓(𝑥) = 𝑥 3 + 3.6𝑥 2 − 36.4

Donde x es la concentración del ion hidrógeno. Encuentre la concentración del ion hidrogeno
para una solución saturada (la acidez es igual a cero).

2. Un nuevo centro de diversiones cuesta $10 millones de pesos y produce una ganancia de $2
millones. Si la deuda se debe pagar en 10 años ¿a qué tasa de interés debe hacerse el préstamo?
El costo actual (P), el pago anual (A) y la tasa de interés (x) se relacionan entre sí mediante la
siguiente formula:

𝑃 (1 + 𝑥)𝑛 − 1
=
𝐴 𝑥 ∗ (1 + 𝑥)𝑛

Sustituyendo datos y simplificando resulta lo siguiente:


(1 + 𝑥)10 − 1
𝑓(𝑥) = −5 = 0
𝑥 ∗ (1 + 𝑥)10

Calcúlese el interés x

3. Suponga que un objeto de masa m se deja caer desde una altura S0 que es la altura del objeto con
respecto al suelo, a los t segundos viene dada por:

𝑚𝑔 𝑚2 𝑔 𝑘𝑡
𝑆(𝑡) = 𝑆0 + 𝑡 − 2 (1 − 𝑒 − 𝑚 )
𝑘 𝑘

Donde S0=300 pies, m=0.25 Slugs, g=-32.17 pies/seg2 y k=0.1 lb x seg/pies.

Obtenga con una precisión dentro de 0.01 seg que es el tiempo que tarda ese objeto en llegar al
suelo.

26
Métodos numéricos

4. Un abrevadero de longitud L tiene una sección transversal en forma de semicírculo con radio r,
cuando se llena de agua hasta una distancia h de la parte superior, el volumen V de agua es:
(Burden y Faires, 2002, p. 55)

𝑉 = 𝐿 [0.5𝜋𝑟 2 − 𝑟 2 𝑎𝑟𝑐𝑠𝑒𝑛 ( ) − ℎ(𝑟 2 − ℎ2 )1/2 ]
𝑟

Suponga que L=10 pies, r=1 pie y que V=12.4 pies3. Determine la profundidad del agua en el
abrevadero hasta 0.01 pies.

5. Realizar ejercicios: 5.13, 5.15 (Chapra y Canale, 2006, p.140)

27
Métodos numéricos

Método de la regla falsa o de la falsa posición


Este método aproxima en forma más eficiente a la raíz. Aprovecha la idea de unir los 2
puntos de intervalo mediante una línea recta en lugar de una curva como en el caso de la
bisección, la intersección de esta línea con el eje x produce una mejor estimación de la raíz,
el remplazamiento de la curva por una recta da una posición falsa de la raíz derivándose de
esto su nombre.

f(x)
f(xb)

xr

(xa)
x
raíz (xb)

f(xa)

El algoritmo es el mismo que el del método de la bisección, lo único que cambia es la


fórmula de xr.

𝒇(𝒙𝒃 ) ∗ (𝒙𝒂 − 𝒙𝒃 )
𝒙𝒓 = 𝒙𝒃 − [𝑅𝑒𝑔𝑙𝑎 𝑓𝑎𝑙𝑠𝑎]
𝒇(𝒙𝒂 ) − 𝒇(𝒙𝒃 )

Ejemplo:
Determinar por el método de la falsa posición la raíz de la siguiente función:
f(x)= e-x –x, tomando 5 decimales.

Es=(0.5 x 102–n)%
Es=(0.5 x 102-5)% = (0.5 x 10-3)% = 0.0005%= 0.05

28
Métodos numéricos

Método por sustitución


x f(x) resultado
-1 3.7183 exp(1) + 1=3.7183
-0.5 2.1487 exp(0.5) + 0.5=2.1487
0 1 exp(0) - 0=1
0.5 0.1065 exp(-0.5) -0.5=0.1065
1 -0.6321 exp(-1) -1=-0.6321
1.5 -1.2769 exp(-1.5) - 1.5=-1.2769

f(xa).f(xr)< f(xa).f(xr)>
f(xa).f(xr)=
iter xa xb f(xa) f(xb) xr f(xr) 0 0 Ea(%) Es(%)
0
xb=xr xa=xr
1 0.5 1.0000 0.1065 -0.6321 0.5721 -0.007779 si no no 100 0.05
2 0.5 0.5721 0.1065 -0.0078 0.5672 -0.000095 si no no 0.865189
3 0.5 0.5672 0.1065 -0.0001 0.5671 -0.000001 si no no 0.010612 parar

Iteración 1
f(xa)=e-0.5 -0.5 = 0.1065
f(xb)=e-1 -1 = -0.6321
xr= 1 - (-0.6321)*(0.5-1)/(0.1065-(-0.6321)) =0.5721
f(xr)=e-0.5721 – 0.5721 =-0.0078 si f(xa)*f(xr) <0 entonces xb=xr
Ea =100%

Iteración 2
f(xa)=e-0.5 -0.5 = 0.1065
f(xb)=e-0.5721 -0.5721 = -0.0078
xr= 0.5721 - (-0.0078)*(0.5-0.5721)/(0.1065-(-0.0078)) =0.5672
f(xr)=e-0.5672 – 0.5672 =-0.00008887 si f(xa)*f(xr) <0 entonces xb=xr
Ea =|(0.5672 - 0.5721)/0.5672| % = 0.0086 % =0.86

Iteración 3
f(xa)=e-0.5 -0.5 = 0.1065
f(xb)=e-0.5672 -0.5672 = -0.000088871
xr= 0.5672 - (-0.000088871)*(0.5-0.5672)/(0.1065-(-0.000088871)) = 0.5671
f(xr)=e-0.5671 – 0.5671 =0.000067843 si f(xa)*f(xr) >0 entonces xa=xr
Ea =|(0.5671 - 0.5672)/0.5671| % = 0.00017634 % =0.01

Algoritmo para el método de la regla falsa o falsa posición


function falsaposicion=falsaposicion()
clc;
f =input('Dame la función : ','s');
xa =input('Dame el intervalo inferior :');
xb =input('Dame el intervalo superior :');
err =input('Dame el porciento de error :');
f =inline(f);

if f(xa)*f(xb)<0
ezplot(f), grid on
erro=100;
xr =0;
i =0;
while erro>err

29
Métodos numéricos

ea=xr;
xr=xb-((f(xb)*(xa-xb))/(f(xa)-f(xb)));
if f(xa)*f(xr)>0
xa=xr;
else
xb=xr;
end
erro=abs(((ea-xr)/xr)*100);
i=i+1;
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',xr,i);
else
fprintf('\n\nNo existe la raíz en el intervalo\n')
end
end

Solución en Matlab:
Dame la función: exp(-x) -x
Dame el intervalo inferior: 0.5
Dame el intervalo superior: 1
Dame el porciento de error: 0.05

Resultado de la raíz: 0.567144 en 3 iteraciones


exp(-x)-x

200

150

100

50

-6 -4 -2 0 2 4 6
x

30
Métodos numéricos

Ejercicios
1. Determinar por el método de la falsa posición la raíz de las siguientes funciones:
1) 𝑓(𝑥) = 6𝑥 3 − 5𝑥 2 + 7𝑥 − 2 tomando 5 decimales
2) 𝑓(𝑥) = 0.7𝑥 5 − 8𝑥 4 + 44𝑥 3 − 90𝑥 2 + 82𝑥 − 25 tomando 5 decimales
3) 𝑓(𝑥) = 𝑠𝑒𝑛 𝑥 – 𝑒 −𝑥 tomando 5 decimales
4) 𝑓(𝑥) = 4 𝑠𝑒𝑛 2𝜋𝑥 + 𝑥 tomando 5 decimales
5) 𝑓(𝑥) = 𝑙𝑜𝑔 (2 + 𝑥) – 𝑥 tomando 5 decimales
6) 𝑓(𝑥) = 0.6 𝑒 −0.3 − 𝑠𝑒𝑛 2𝑥 tomando 5 decimales
7) 𝑓(𝑥) = 𝑥 𝑙𝑜𝑔 𝑥 − 10 en el intervalo [5,6] tomando 5 decimales
8) 𝑓(𝑥) = 2𝑥 0.2 − 𝑒 −𝜋𝑥 𝑡𝑎𝑛(𝑥) − 2 en el intervalo [1, 1.55] tomando 5 decimales
9) 𝑓(𝑥) = 2𝑥 0.6 − 𝑐𝑜𝑠(𝑥) 𝑙𝑜𝑔10 (𝑥) − 20 en el intervalo [40,55] tomando 5 decimales
10) 𝑓(𝑥) = 𝑥 3 𝑒 −𝑥 + 4𝑥 2 − 10 en el intervalo [1, 2] tomando 5 decimales

2. Determinar el factor de fricción f para los flujos turbulentos en una tubería que está dada por:
1 𝑒 9.35
= 1.14 − 2𝑙𝑜𝑔10 ( + )
√𝑓 𝐷 𝑅𝑒√𝑓

Llamada correlación de Colebrook, donde Re es el número de Reynolds, e es la aspereza de la


superficie de la tubería y D es el diámetro de la tubería.
Tomemos D =0.1 m, e =0.0025 m y Re =30000.

3. Para el flujo turbulento de un fluido a través de un tubo liso, es posible establecer la siguiente
relación entre el factor de fricción cf y el número de Reynolds Re:
1
√ = −0.4 + 1.74 ln(𝑅𝑒 √𝑐𝑓)
𝑐𝑓
Calcular cf para Re=104, 105 y 106.

4. Basado en el trabajo de Frank-Kamenetski en 1955, las temperaturas en el interior de un material


con fuentes de calor incrustadas pueden determinarse si resolvemos esta ecuación: (Curtis y
Wheatly, 2004, p.72).
1
−( )𝑡
1
( )𝑡 1
𝑒 2 𝑐𝑜𝑠ℎ−1 (𝑒 2 ) = √ 𝐿𝑐𝑟
2
Dado que 𝐿𝑐𝑟 = 0.088, encuentre 𝑡.
3
5. Una esfera de densidad 𝑑 y radio 𝑟 pesa 4
𝜋𝑟 3 𝑑. El volumen de un segmento esférico es
1
3
𝜋(3𝑟ℎ2 − ℎ3 ).

Encuentre la profundidad a la que una esfera de densidad 0.6 se hunde en el agua como una
fracción de su radio (Ver la figura adjunta) (Curtis y Wheatly, 2004, p.74).

31
Métodos numéricos

Métodos abiertos
Se basan en fórmulas que requieren únicamente de un solo valor de inicio x o que empiecen
con un par de ellos, pero que no necesariamente encierran a la raíz. Como tales, algunas
veces divergen o se aleja de la raíz verdadera a medida que crece el número de iteraciones.
Sin embargo, cuando estos métodos convergen por lo general lo hacen mucho más rápido
que los métodos que usan intervalos.

Método del punto fijo


Se emplea una fórmula para predecir la raíz. Esta fórmula puede desarrollarse como una
iteración simple de punto fijo, al arreglar la ecuación 𝑓(𝑥) = 0 de tal modo que x esté del
lado izquierdo de la ecuación.

𝑥 = 𝑔(𝑥)

Esta transformación se puede llevar a cabo mediante operaciones algebraicas o


simplemente agregando x a cada lado de la ecuación original.

Ejemplo:
𝑥2+ 3
𝑥 2 − 2𝑥 + 3 = 0 Se arregla para obtener 𝑥 = (operación algebraica)
2

Ejemplo:
𝑠𝑒𝑛 𝑥 = 0 Se arregla para obtener 𝑥 = 𝑠𝑒𝑛 𝑥 + 𝑥 (sumando x en ambos lados)

La ecuación x=g(x) nos proporciona una fórmula para predecir un nuevo valor de x en
función del valor anterior de x, se utiliza para obtener una nueva aproximación xi+1
expresada por la formula iterativa:

𝒙𝒊+𝟏 = 𝒈(𝒙𝒊 ) [𝑃𝑢𝑛𝑡𝑜 𝑓𝑖𝑗𝑜]

32
Métodos numéricos

El error aproximado porcentual está dado por:


𝑥𝑖+1 − 𝑥𝑖
𝐸𝑎 = | | 100% [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑥𝑖+1

y y=x
g(x1)

x2=g(x1)

g(x2)
x3=g(x2)
y=g(x)
x
x1 x3 x5 x6 x4 x2

Ejemplo:
Determinar por el método del punto fijo la raíz de la siguiente función:
𝑓(𝑥) = 𝑒 −𝑥 − 𝑥, con punto de inicio en 0 y tomando 5 decimales.

𝑒 −𝑥 − 𝑥 = 0
Se arregla para obtener 𝑥 = 𝑒 −𝑥
iter xi xi+1 Ea(%) Es(%)
1 0 1 100 0.05
2 1 0.36787944 171.828183
3 0.36787944 0.69220063 46.8536395
4 0.69220063 0.5004735 38.3091466
5 0.5004735 0.60624354 17.4467897
6 0.60624354 0.54539579 11.1566225
7 0.54539579 0.57961234 5.90335081
8 0.57961234 0.56011546 3.48086698
9 0.56011546 0.57114312 1.93080393

33
Métodos numéricos

10 0.57114312 0.56487935 1.10886824


11 0.56487935 0.56842873 0.62441912
12 0.56842873 0.56641473 0.35556841
13 0.56641473 0.56755664 0.20119652
14 0.56755664 0.56690891 0.11425564
15 0.56690891 0.56727623 0.06475157
16 0.56727623 0.5670679 0.03673877 Parar

Nota: Cada iteración se acerca cada vez más al valor estimado con el valor verdadero de la
raíz a 0.5670679, continuando las iteraciones hasta que Ea<= 0.05

Iteración 1
f(x(0))=e-0 = 1
Ea=100

Iteración 2
f(x(1))= e-1=0.3679
Ea=|(0.3679-1)/(0.3679)|%= 1.7219% =172.19

Iteración 3
f(x(0.3679))= e-0.3679 =0.6922
Ea=|(0.6922-0.3679)/(0.6922)|%=0.4685%=46.85

Iteración 4
f(x(0.6922))= e-0.6922=0.5005
Ea=|(0.5005-0.6922)/(0.5005)|% =0.3830% = 38.30

Iteración 5
f(x(0.5005))= e-0.5005=0.6062
Ea=|(0.6062 -0.5005)/(0.6062)|% = 0.1744% = 17.44

Iteración 6
f(x(0.6062))= e-0.6062 = 0.5454
Ea=|(0.5454 - 0.6062)/(0.5454)|% = 0.1115% = 11.15

Iteración 7
f(x(0.5454))= e-0.5454= 0.5796
Ea=|(0.5796-0.5454)/(0.5796)|% = 0.0590% = 5.90

Iteración 8
f(x(0.5796))= e-0.5796= 0.5601
Ea=|(0.5601 – 0.5796)/(0.5601)|% = 0.0348% = 3.48
Iteración 9
f(x(0.5601))= e-0.5601 = 0.5712
Ea=|(0.5712 – 0.5601)/(0.5712)|% =0.0194% = 1.94

Iteración 10
f(x(0.5712))= e-0.5712 =0.5648
Ea=|(0.5648 – 0.5712)/(0.5648)|% =0.0113% = 1.13

34
Métodos numéricos

Iteración 11
f(x(0.5648)) = e-0.5648 = 0.5685
Ea=|(0.5685 -0.5648)/(0.5685)|% = 0.0065% = 0.65

Algoritmo para el método del punto fijo


function puntofijo=puntofijo()
clc;
f =input('Dame la función : ','s');
xi =input('Dame el punto de inicio :');
err=input('Dame el porciento de error :');
f =inline(f);

ezplot(f), grid on
i =0;
ea =100;
while err <= ea
xr=f(xi);
ea=abs(((xr-xi)/xr)*100);
xi=xr;
i=i+1;
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',xr,i);
end

Solución en Matlab:
Dame la función: exp(-x)
Dame el punto de inicio: 1
Dame el porciento de error: 0.05

Resultado de la raíz= 0.567068 en 15 iteraciones

35
Métodos numéricos

Ejercicios
1. Determinar por el método del punto fijo la raíz de las siguientes funciones:
1) 𝑓(𝑥) = 𝑠𝑒𝑛 𝑥, con punto de inicio en 1 y tomando 5 decimales
2) 𝑓(𝑥) = 𝑐𝑜𝑠 𝑥 − 3𝑥, con punto de inicio en 0 y tomando 5 decimales
3) 𝑓(𝑥) = 𝑥 4 + 8𝑥 3 + 11𝑥 2 − 32𝑥 − 60, con punto de inicio en 1 y tomando 5
decimales
4) 𝑓(𝑥) = 8𝑥 + 𝑐𝑜𝑠(𝑥 + 𝜋) − 𝑥 2 con punto de inicio en 2 y tomando 5 decimales
5) 𝑓(𝑥) = 𝑥 3 𝑐𝑜𝑠(𝑥) − 5𝑥 2 − 1 en el intervalo [36, 37] y tomando 5 decimales
6) 𝑓(𝑥) = 1000𝑥 2 𝑡𝑎𝑛(𝑥) 𝑒 −5𝑥 − 10 en el intervalo [0.5, 1] y tomando 5 decimales
7) 𝑓(𝑥) = 2𝑥 − 20𝑠𝑒𝑛(𝑥) iniciando en el intervalo [2,4] y tomando 5 decimales
8) 𝑓(𝑥) = 𝑥 2 𝑒 −𝑥 + 𝑙𝑛𝑥 – 3 con punto de inicio en 1 y tomando 5 decimales
9) 𝑓(𝑥) = 𝑠𝑒𝑛(0.1𝑥) + 𝑒 −2𝑥 − 0.6, en el intervalo [20, 30] y tomando 5 decimales
10) 𝑓(𝑥) = 𝑐𝑜𝑠(3𝑥) + 5𝑒 −0.01𝑥 − 3, en el intervalo [34, 34.5] y tomando 5 decimales.

2. En estudios de recolección de energía solar al enfocar un campo de espejos planos en un


colector central, un investigador obtuvo esta ecuación para el factor de concentración
geométrico 𝐶: (Curtis y Wheatly, 2004, p.72).
ℎ 2
𝜋( ) 𝐹
cos 𝐴
𝐶=
0.5𝜋𝐷2 (1 + 𝑠𝑒𝑛𝐴 − 0.5𝑐𝑜𝑠𝐴)

donde 𝐴 es el ángulo del borde del campo, 𝐹 es la cobertura fraccional del campo con
espejos, 𝐷 es el diámetro del colector y ℎ es la altura del colector. Encuentre 𝐴 si ℎ = 300,
𝐶 = 1200, 𝐹 = 0.8, y 𝐷 = 14.

3. DeSantis (1976) ha derivado una relación para el factor de compresibilidad de gases reales
de la forma: (Curtis y Wheatly, 2004, p.72).

1 + 𝑦 + 𝑦2 − 𝑦3
𝑧=
(1 − 𝑦)3

donde 𝑦 = 𝑏/(4𝑣), siendo b la corrección de Van der Waals y 𝑣 el volumen molar. Si


𝑧 = 0.892, ¿cuál es el valor de 𝑦?

4. Lee y Duffy (1976) relacionan el factor de fricción para el flujo de una suspensión de
partículas fibrosas con el número de Reynolds mediante esta ecuación empírica: (Curtis y
Wheatly, 2004, p.72).
1 1 5.6
= ( ) 𝐼𝑛(𝑅𝐸√𝑓) + (14 − )
√𝑓 𝑘 𝑘
En su relación, 𝑓 es el factor de fricción, 𝑅𝐸 es el número de Reynolds y 𝑘 es una constante
determinada por la concentración de la suspensión. Para una suspensión con una
concentración de 0.08%, 𝑘 = 0.28. ¿Cuál es el valor de 𝑓 si 𝑅𝐸 = 3750?

5. Realizar ejercicio: 23 (Burden y Faire, 2002, p.65)

36
Métodos numéricos

Método de Newton-Raphson
Este método es el más utilizado. Si el valor inicial para la raíz es xi, entonces se puede
trazar una tangente desde el punto [xi, f(xi)] de la curva, el punto donde esta tangente cruza
al eje x representa una aproximación mejorada de la raíz.

f(x)
Pendiente=f’(xi)

f(xi)

f(xi)-0

x
xi+1 xi

xi -xi+1

El método de Newton-Raphson se deduce a partir de esta interpretación geométrica y se


tiene que la primera derivada en x es equivalente a la pendiente:
f(xi ) − 0
f ′ (xi ) =
xi − xi+1

Que se arregla para obtener

𝐟(𝐱𝐢 )
𝐱𝐢+𝟏 = 𝐱𝐢 – [𝑁𝑒𝑤𝑡𝑜𝑛 − 𝑅𝑎𝑝ℎ𝑠𝑜𝑛]
𝐟 ′ (𝐱𝐢 )

El error aproximado porcentual está dado por:

𝑥𝑖+1 − 𝑥𝑖
𝐸𝑎 = | | 100% [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑥𝑖+1

37
Métodos numéricos

Ejemplo:
Determinar por el método del Newton-Raphson la raíz de la siguiente función:

𝑓(𝑥) = 𝑒 −𝑥 − 𝑥, con punto de inicio en 0 y tomando 5 decimales.

Derivada de la función: 𝑓’(𝑥) = – 𝑒 −𝑥 − 1

iter xi f(xi) f'(xi) xi+1 Ea(%) Es(%)


1 0 1 -2 0.5 100 0.05
2 0.5 0.10653066 -1.60653066 0.566311 11.709291
3 0.566311 0.00130451 -1.56761551 0.56714317 0.14672871
4 0.56714317 1.9648E-07 -1.56714336 0.56714329 0.000022 parar

Iteración 1
xi=0
f(xi)= e0 -0= 1
f’(xi)=-e0-1=-2
xi+1=0-(1/-2)= 0.50
Ea=|(0.50 - 0)/(0.50)|% =1 % =100

Iteración 2
xi=0.50
f(xi)= e-0.50 -0.50= 0.1065
f’(xi)=-e-0.50-1=-1.6065
xi+1=0.50-(0.1065/-1.6065)= 0.5663
Ea=|(0.5663–0.50)/(0.5663)|% =0.1171 % =11.71

Iteración 3
xi=0.5663
f(xi)= e-0.5663 -0.5663= 0.0013
f’(xi)=-e-0.5663-1=-1.5676
xi+1=0.5663-(0.0013/-1.5676)= 0.5671
Ea=|(0.5671 – 0.5663)/(0.5671)|% =0.0014 % =0.14

Iteración 4
xi=0.5671
f(xi)= e-0.5671 -0.5671=0.0000678
f’(xi)=-e-0.5671-1=-1.5672
xi+1=0.5671-(0.0000678/-1.5672)= 0.5671
Ea=|(0.5671 – 0.5671)/(0.5671)|% =0.00% =0.00

Algoritmo para el método de Newton-Raphson


function newtonraphson=newtonraphson()
clc;
% syms x
f =input('Dame la función : ','s');
dx =input('Dame la derivada de la función: ','s');
%dx =diff(f);

38
Métodos numéricos

pi =input('Dame el valor del punto de inicio: ');


err =input('Dame el porciento de error :');
f =inline(f);
ezplot(f), grid on
dx =inline(dx);
ea =100;
i =0;
while ea>err
xi=pi-(f(pi)/dx(pi));
ea=abs(((xi-pi)/xi)*100);
pi=xi;
i =i+1;
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',pi,i);
end

Solución en Matlab:
Dame la función: 𝑒𝑥𝑝(−𝑥) − 𝑥
Dame la derivada de la función: −𝑒𝑥𝑝(−𝑥) − 1
Dame el valor del punto de inicio: 0
Dame el porciento de error: 0.05

Resultado de la raíz= 0.567143 en 4 iteraciones

39
Métodos numéricos

Método de Newton-Raphson modificado para el cálculo de raíces múltiples


El método modificado utiliza una segunda derivada y la siguiente ecuación.

𝒇(𝒙𝒊 ) ∗ 𝒇′ (𝒙𝒊 )
𝒙𝒊+𝟏 = 𝒙𝒊 − [𝑁𝑒𝑤𝑡𝑜𝑛 − 𝑅𝑎𝑝ℎ𝑠𝑜𝑛 𝑚𝑜𝑑𝑖𝑓𝑖𝑐𝑎𝑑𝑜]
[𝒇′(𝒙𝒊 )]𝟐 − (𝒇(𝒙𝒊 ) ∗ 𝒇′′ (𝒙𝒊 ))

𝑓’’(𝑥𝑖) = segunda derivada de la función

Ejemplo:
Determinar por el método del Newton-Raphson Modificado la raíz de la siguiente función
𝑓(𝑥) = 𝑒 −𝑥 − 𝑥 con punto de inicio en 0 y tomando 5 decimales.

1ª Derivada de la función: 𝑓’(𝑥) = −𝑒 −𝑥 − 1


2ª Derivada de la función 𝑓’’(𝑥) = 𝑒 −𝑥

iter xi f(xi) f'(xi) f''(xi) xi+1 Ea(%) Es(%)


1 0 1.000000 -2.000000 1.000000 0.666667 100 0.05
2 0.666667 -0.153250 -1.513417 0.513417 0.568769 17.212195
3 0.568769 -0.002547 -1.566222 0.566222 0.567144 0.286570
4 0.567144 -0.000001 -1.567143 0.567143 0.567143 0.000084 parar

Algoritmo para el método de Newton-Raphson Modificado


function newtonraphsonmodificado=newtonraphsonmodificado()
clc;
f =input('Dame la función : ','s');
dx =input('Dame la 1a. derivada de la función: ','s');
dx2 =input('Dame la 2a. derivada de la función: ','s');
pi =input('Dame el valor del punto de inicio: ');
err =input('Dame el porciento de error :');
f =inline(f);

ezplot(f), grid on
dx =inline(dx);
dx2 =inline(dx2);
ea =100;
i =0;
while ea>err
xi=pi-(f(pi)*dx(pi))/((dx(pi)^2)-(f(pi)*dx2(pi)));
ea=abs(((xi-pi)/xi)*100);
pi=xi;
i =i+1;
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',pi,i);
end

40
Métodos numéricos

Solución en Matlab:
Dame la función f(x) : 𝑒𝑥𝑝(−𝑥) − 𝑥
Dame la 1a. derivada de la función f(x) : −𝑒𝑥𝑝(−𝑥) − 1
Dame la 2a. derivada de la función f(x) : 𝑒𝑥𝑝(−𝑥)
Dame el valor inicial de x: 0
Dame el porciento del error : 0.05

Resultado de la raíz= 0.567143 en 4 iteraciones

41
Métodos numéricos

Ejercicios
1. Determinar por el método del Newton Raphson y el método de Newton Raphson
modificado la raíz de las siguientes funciones:

1. f(x) = 6x3 -8x2 -10x+3 con punto de inicio en 0 y tomando 5 decimales


2. f(x) = x2 – 3x + 2 - eX con punto de inicio en 1 y tomando 5 decimales
3. f(x) = x2 -2xe-x +e-2x en el intervalo [0,1] y tomando 5 decimales
4. f(x) = x2cosx − 6x lnx − 25 con punto de inicio en 23 y tomando 5 decimales
5. f(x) = (x − 2)2 – ln x con punto de inicio en 1 y tomando 5 decimales
6. f(x) = cos(x)− 3x con punto de inicio en 0.5 y tomando 5 decimales
7. f(x) = cos(5x) + 5sen(15x) +e−0.05x − 3 en el intervalo [16.1, 16.5] y tomando 5
decimales
8. f(x) = log10(50 x − 40) – x1.2 + 20 cos(x) + 12 en el intervalo [5, 10] y tomando 5
decimales

2. Resuelva 𝑥 3 − 2𝑥 − 5 = 0. Esta ecuación tiene valor histórico: fue la ecuación que usó
John Wallis para presentar por primera vez el método de Newton a la academia francesa de
ciencias en el siglo XV.

3. La tasa de flujo de agua a través de una corriente a menudo se mide instalando un


vertedero. Esto equivale a construir una presa a través del arroyo con una muesca en forma
de V cerca del centro (el punto de la muesca es abajo). Si se descuida la velocidad aguas
𝑓𝑡 3
arriba, el flujo 𝑄 (𝑠𝑒𝑔) está relacionado con la distancia ℎ (pies) desde la superficie del agua
arriba hasta el punto de la V y con el ángulo ϴ (grados) entre los lados de la muesca por
esta fórmula: (Curtis y Wheatly, 2004, p.75).
8 𝛳
𝑄 = 0.59 ∗ ( ) ∗ 𝑡𝑎𝑛 ( ) ∗ √(2𝑔) ∗ ℎ2.5
15 2
𝑓𝑡
donde 𝑔 es la constante gravitacional 32.2 𝑠𝑒𝑐 2 .Si 𝑄 = 200, haga una tabla que muestre
cómo ℎ se relaciona con ϴ para valores de ϴ entre 20 y 130 grados.

4. La velocidad v de un cohete Saturno en vuelo vertical cerca de la superficie de la Tierra


puede ser aproximada por:

𝑀0
𝑣 = 𝑢𝐿𝑛 − 𝑔𝑡
𝑀0 − 𝑚𝑡

Donde:
U = 2510 m/s =velocidad de escape relativa al cohete
M0 = 2.8 x106 kg = masa del cohete al despegar
m = 13.3 x 103 kg/s = tasa de consumo de combustible
g = 9.81 m/s2 = aceleración gravitacional
t = tiempo medido desde el despegue
Determinar el tiempo en que el cohete alcanza la velocidad del sonido (335 m/s)

42
Métodos numéricos

Método de la Secante
En el método de Newton-Raphson el problema que existe es la evaluación de la derivada
por lo que en este método en lugar de una derivada se utiliza una diferencia dividida finita
regresiva.

f(x)

f(xi)

f(xi-1)

x
xi+1 xi

El planteamiento requiere de dos puntos iniciales de x, sin embargo, debido a que no


requiere que f(x) cambie de signo entre estos valores, este método no es clasificado como
aquellos que usan intervalos.

𝒇(𝒙𝒊 ) ∗ (𝒙𝒊−𝟏 − 𝒙𝒊 )
𝒙𝒊+𝟏 = 𝒙𝒊 − [𝑆𝑒𝑐𝑎𝑛𝑡𝑒]
𝒇(𝒙𝒊−𝟏 ) − 𝒇(𝒙𝒊 )
El error aproximado porcentual está dado por:

𝑥𝑖+1 − 𝑥𝑖
𝐸𝑎 = | | 100% [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑥𝑖+1

Ejemplo:
Determinar por el método de la Secante la raíz de la siguiente función:
𝑓(𝑥) = 𝑒 −𝑥 − 𝑥, con valores iniciales en xi-1=0 y xi=1 y tomando 5 decimales.

43
Métodos numéricos

iter xi-1 xi f(xi-1) f(xi) xi+1 Ea(%) Es(%)


1 0 1 1.000000 -0.632121 0.612700 100 0.05
2 1.000000 0.612700 -0.632121 -0.070814 0.563838 8.665860
3 0.612700 0.563838 -0.070814 0.005182 0.567170 0.587472
4 0.563838 0.567170 0.005182 -0.000042 0.567143 0.004770 parar

Iteración 1
xi-1=0
xi=1
f(xi-1)= e0 -0= 1
f(xi)=e-1-1= -0.6321
xi+1=1-((-0.6321)(0-1))/(1-(-0.6321))=0.6127
Ea=|(0.6127 - 1)/(0.6127)|% =0.6127 % =61.27

Iteración 2
xi-1=1
xi=0.6127
f(xi-1)= e-1 -1= -0.6321
f(xi)=e0.6127-0.6127=-0.0708
xi+1=0.6127-((-0.0708)(1-0.6127)/(-0.6321-(-0.0708))=0.5638
Ea=|(0.5638-0.6127)/(0.5638)|% =0.0867 % =8.67

Iteración 3
xi-1=0.6127
xi=0.5638
f(xi-1)= e-0.6127-0.6127=-0.0708
f(xi)=e-0.5638-0.5638=0.0052
xi+1=0.5638-((0.0052)(0.6127-0.5638)/(-0.0708-(0.0052))=0.5671
Ea=|(0.5671-0.5638)/(0.5671)|% =0.0058 % =0.58

Iteración 4
xi-1=0.5638
xi=0.5671
f(xi-1)= e-0.5638-0.5638=0.0052
f(xi)=e-0.5671-0.5671=0.0000678
xi+1=0.5671-((0.0000678)(0.5638-0.5671)/(-0.0052-(0.0000678))=0.5671
Ea=|(0.5671-0.5671)/(0.5671)|% = 0.00% =0.0004

Algoritmo para el método de la Secante


function secante=secante()
clc;
f =input('Dame la función : ','s');
pa =input('Dame el punto para xi-1:');
pb =input('Dame el punto para xi :');
err =input('Dame el porciento de error :');
f =inline(f);

ezplot(f), grid on
ea =100;

44
Métodos numéricos

i =0;
while ea>err
xi =pb-((f(pb)*(pa-pb))/(f(pa)-f(pb)));
ea =abs(((xi-pb)/xi)*100);
pa =pb;
pb =xi;
i =i+1;
end
fprintf('\n\nResultado de la raíz=%12.6f en %4d iteraciones\n',xi,i);
end

Solución en Matlab:
Dame la función: exp(-x)-x
Dame el punto xi-1: 0
Dame el punto xi : 1
Dame el porciento de error: 0.05
Resultado de la raíz= 0.567143 en 4 iteraciones

45
Métodos numéricos

Ejercicios
1. Determinar por el método de la secante la raíz de las siguientes funciones:

1) f(x) = x3 + 2x2 + 10x - 20 con valores inicial en xi-1=0 y xi=1 y tomando 5 decimales
2) f(x) = x3 - 2x2cos(2x) - 12 con valores inicial en xi-1=0 y xi=4 y tomando 5 decimales
3) f(x) = x4 + 8x3 + 11x2 - 32x - 60 con valor inicial en xi-1=1 y xi=3 y tomando 5 decimales
4) f(x) = 3xsen(x) – ex con valor inicial en xi-1=1.5 y xi=0.8 y tomando 5 decimales
5) 𝑓(𝑥) = 𝑥 – 2𝑐𝑜𝑠(𝑥) con valor inicial en xi-1 =1 y xi =1.5 y tomando 5 decimales
6) f(x) = 5xe−x + cos(5x) con valor inicial en xi-1= 3.9 y xi = 4 y tomando 5 decimales
7) f(x) = x2 cos(10x) + 1 con valor inicial en xi-1 =1.3 y xi =1.4 y tomando 5 decimales
8) f(x) = 5x cos(x) + 2 con valor inicial en xi-1 =1.0 y xi =1.1 y tomando 5 decimales
9) f(x) = e−2x sen (0.1x) − e−10x cos (0.1x) + 0.3 con valor inicial en xi-1 = 0.4 y xi = 0.3 y
tomando 5 decimales
10) f(x)= 2000x3e-5xcos(2x) + 1 con valor inicial en xi-1= 1.5 y xi=1.7 tomando 5 decimales

2. La siguiente ecuación permite calcular la concentración de un químico en un reactor donde se


tiene una mezcla completa:

𝑐 = 𝑐𝑒𝑛𝑡 (1 – 𝑒 –0.04𝑡 ) + 𝑐0 𝑒 –0.04𝑡

Si la concentración inicial es 𝑐0 = 5 y la concentración de entrada es 𝑐𝑒𝑛𝑡 = 12, calcule el


tiempo requerido para que 𝑐 sea el 85% de 𝑐𝑒𝑛𝑡 .

3. El valor acumulado de una cuenta de ahorros que se basa en pagos periódicos puede calcularse
con la ecuación de anualidad vencida.

𝑃
𝐴= [(1 + 𝑖)𝑛 − 1]
𝑖

En esta ecuación, A es el monto de la cuenta, P es la cantidad que se deposita periódicamente e


i es la tasa de interés por periodo para los n periodos de depósito. A un ingeniero le gustaría
tener una cuenta de ahorros con un monto de 750,000 dólares al momento de retirarse dentro
de 20 años, y puede depositar 1,500 dólares mensuales para lograr dicho objetivo. ¿Cuál es la
tasa mínima de interés a que puede invertirse ese dinero, suponiendo que es un interés
compuesto mensual?

4. Realizar ejercicio: 17,18 (Kiusalaas, 2010, p.166)

46
Métodos numéricos

Método de Müller

Consiste en obtener los coeficientes de la parábola que pasa por tres puntos elegidos.
Cuyos coeficientes son sustituidos en la formula cuadrática para obtener el valor donde la
parábola intercepta al eje x; es decir, la raíz estimada. La aproximación se puede facilitar, si
se escribe la ecuación de la parábola en una forma conveniente. Una de las mayores
ventajas del método de Müller, es que al trabajar con la formula cuadrática es posible
encontrar las raíces reales, tanto como las raíces complejas. El método de Müller en si es
una generalización del método de la Secante.

p(x)=a0 + a1x + a2x2

y Parábola
y=f(x)

Raíz

x
x2 x1 x0

Raíz estimada

En el método de Müller se usan tres aproximaciones iníciales x0, x1 y x2 con las cuales
procederíamos a determina la siguiente aproximación x3, considerando la intercepción del
eje x con la parábola que pasa por (x0, f(x0)), (x1, f(x1) y (x2, f(x2)).

−𝟐𝒄
𝒙𝟑 = 𝒙𝟐 + [𝑀𝑢𝑙𝑙𝑒𝑟]
𝒃 ± √𝒃𝟐 − 𝟒𝒂𝒄

Dónde:

47
Métodos numéricos

𝒄 = 𝑓(𝑥2 )

(𝑥0 − 𝑥2 )2 (𝑓(𝑥1 ) − 𝑓(𝑥2 )) − (𝑥1 − 𝑥2 )2 (𝑓(𝑥0 ) − 𝑓(𝑥2 ))


𝒃=
(𝑥0 − 𝑥2 )(𝑥1 − 𝑥2 )(𝑥0 − 𝑥1 )

(𝑥1 − 𝑥2 )(𝑓(𝑥0 ) − 𝑓(𝑥2 )) − (𝑥0 − 𝑥2 )(𝑓(𝑥1 ) − 𝑓(𝑥2 ))


𝒂=
(𝑥0 − 𝑥2 )(𝑥1 − 𝑥2 )(𝑥0 − 𝑥1 )

El error aproximado porcentual está dado por:

𝑥3 − 𝑥2
𝐸𝑎 = | | 100% [𝐸𝑟𝑟𝑜𝑟 𝑎𝑝𝑟𝑜𝑥𝑖𝑚𝑎𝑑𝑜 𝑝𝑜𝑟𝑐𝑒𝑛𝑡𝑢𝑎𝑙]
𝑥3

Condición:

𝑆𝑖 |𝑏 + √𝑏 2 − 4𝑎𝑐| > |𝑏 − √𝑏 2 − 4𝑎𝑐| 𝑒𝑛𝑡𝑜𝑛𝑐𝑒𝑠: 𝑏 + √𝑏 2 − 4𝑎𝑐,


𝑠𝑖 𝑛𝑜: 𝑏 − √𝑏 2 − 4𝑎𝑐

Realizar el método de Müller en Excel y en Matlab.

48
Métodos numéricos

Comparación de las características de los métodos alternativos para encontrar raíces de


ecuaciones algebraicas y transcendentales. Las comparaciones se basan en la experiencia
general y no toman en cuenta el comportamiento de funciones específicas.

Tabla 1. Comparación de los métodos para encontrar raíces (Elaborado a partir de Chapra y Canale,
2006, p.227)

Amplitud
Valores Velocidad de Complejidad de
Método Estabilidad Exactitud de Comentarios
iniciales convergencia programación
aplicación
Analítico - - - - Limitada
Puede
tomar más
Raíces tiempo
Grafico - - - Pobre -
reales que el
método
numérico
Raíces
Bisección 2 Lenta Siempre Buena Fácil
reales
Falsa Raíces
2 Lenta/media Siempre Buena Fácil
posición reales
Posiblemente
Punto fijo 1 Lenta Buena General Fácil
divergente
Requiere la
Newton- Posiblemente
1 Rápida Buena General Fácil evaluación
Raphson divergente
de ƒ′(x)
Rápida para
Requiere la
Newton- raíces
Posiblemente evaluación
Raphson 1 múltiples; Buena General Fácil
divergente de ƒ″(x) y
modificado media
ƒ′(x)
para una sola
Los valores
Posiblemente iniciales
Secante 2 Media a rápida Buena General Fácil
divergente no tiene que
acotar la raíz

Posiblemente
Müller 3 Media a rápida Buena Polinomios Moderada
divergente

49
Métodos numéricos

SOLUCIÓN DE ECUACIONES ALGEBRAICAS LINEALES SIMULTÁNEAS


Se denomina ecuación lineal a aquella que tiene la forma de un polinomio de primer grado,
es decir, las incógnitas no están elevadas a potencias, ni multiplicadas entre sí, ni en el
denominador.
Por ejemplo, 3𝑥 + 2𝑦 + 6𝑧 = 6 es una ecuación lineal con tres incógnitas.
Las ecuaciones lineales con 2 incógnitas representan una recta en el plano y si la ecuación
lineal tiene 3 incógnitas, su representación gráfica es un plano en el espacio como se
muestra en la Figura 4.

Figura 4. Representación gráfica de ecuación lineal con 2 y 3 incógnitas

Solución de un sistema de ecuaciones: es un conjunto de valores de las incógnitas que


verifican simultáneamente a todas y cada una de las ecuaciones del sistema.

Tipos de sistemas de ecuaciones lineales

𝐼𝑛𝑐𝑜𝑚𝑝𝑎𝑡𝑖𝑏𝑙𝑒 (𝑁𝑜 𝑡𝑖𝑒𝑛𝑒 𝑠𝑜𝑙𝑢𝑐𝑖ó𝑛)


{ 𝐷𝑒𝑡𝑒𝑟𝑚𝑖𝑛𝑎𝑑𝑜𝑠 (𝑺𝒐𝒍𝒖𝒄𝒊ó𝒏 ú𝒏𝒊𝒄𝒂)
𝐶𝑜𝑚𝑝𝑎𝑡𝑖𝑏𝑙𝑒 (𝑆𝑖 𝑡𝑖𝑒𝑛𝑒 𝑠𝑜𝑙𝑢𝑐𝑖ó𝑛) {
𝐼𝑛𝑑𝑒𝑡𝑒𝑟𝑚𝑖𝑛𝑎𝑑𝑜𝑠 (𝐼𝑛𝑓𝑖𝑛𝑖𝑡𝑎𝑠 𝑠𝑜𝑙𝑢𝑐𝑖𝑜𝑛𝑒𝑠)

50
Métodos numéricos

Figura 5. Representación gráfica de los tipos de sistemas de ecuaciones lineales (Kolman y Hill,
2006).

Existen diversos métodos alternativos para encontrar soluciones de ecuaciones algebraicas


lineales simultáneas como es el caso de los métodos (el método Grafico y la Regla de
Cramer) que están limitados a pocas ecuaciones de modo que tienen escasa utilidad para
resolver problemas prácticos.
Un sistema de m ecuaciones lineales en n incógnitas tiene la forma general:

𝑎1,1 𝑥1 + 𝑎1,2 𝑥2 + … 𝑎1,𝑛 𝑥𝑛 = 𝑏1


𝑎2,1 𝑥1 + 𝑎2,2 𝑥2 + … 𝑎2,𝑛 𝑥𝑛 = 𝑏2
𝑎3,1 𝑥1 + 𝑎3,2 𝑥2 + … 𝑎3,𝑛 𝑥𝑛 = 𝑏3
⋮ ⋮ ⋱ ⋮ ⋮
𝑎𝑚,1 𝑥1 + 𝑎𝑚,2 𝑥2 + … 𝑎𝑚,𝑛 𝑥𝑛 = 𝑏𝑚

Con la notación matricial la ecuación se puede escribir como Ax=b, donde A es la matriz
de coeficientes, x es el vector de incógnitas y b es el vector de términos del lado derecho.

𝑎1,1 𝑎1,2 … 𝑎1,𝑛 𝑥1 𝑏1


𝑎2,1 𝑎2,2 … 𝑎2,𝑛 𝑥2 𝑏2
𝐴 = 𝑎3,1 𝑎3,2 … 𝑎3,𝑛 𝑥 = 𝑥3 𝑏 = 𝑏3
⋮ ⋮ ⋱ ⋮ ⋮ ⋮
[𝑎𝑚,1 𝑎𝑚,2 … 𝑎𝑚,𝑛 ] [𝑥𝑚 ] [𝑏𝑚 ]

Los métodos numéricos que se estudian para la solución de un sistema de ecuaciones


lineales se clasifican en dos tipos: directos e iterativos.

Los métodos directos nos proporcionan una solución del sistema en un número finito de
pasos; se usa aritmética finita para los cálculos, se obtiene por lo general una solución
aproximada, debido únicamente a los errores de redondeo, puesto que no hay errores de

51
Métodos numéricos

truncamiento o de fórmula. Los métodos directos más usados tienen como base la
eliminación de Gauss.

En los métodos iterativos se parte de una aproximación inicial a la solución del sistema
dado y se genera, a partir de dicha aproximación, una sucesión de vectores que si converge
lo hace a la solución del sistema. Tendremos fórmulas para calcular los términos de la
sucesión, así que en general no se espera calcular el límite de la sucesión, por lo que
debemos tomar algún término de la sucesión como una solución aproximada del sistema.
Esta vez, además de los errores de redondeo si se usa aritmética finita, habrá errores de
truncamiento o de fórmula. Los métodos iterativos más simples y conocidos están basados
en iteraciones de Punto Fijo.

Métodos directos
Método de eliminación de Gauss
Es un proceso que convierte a la matriz de coeficientes A de n x m en una matriz triangular
superior, mediante la aplicación sistemática de transformaciones elementales de renglón.
Una vez obtenida la matriz triangular superior se aplica un procedimiento conocido como
sustitución hacia atrás para obtener el vector solución x.

Las transformaciones elementales de renglón son:


1. La fila i de la matriz puede ser multiplicada por una constante 𝜆 ≠ 0.
𝜆𝑅𝑖 = 𝑅𝑖
2. A la fila i de una matriz le puede ser sumada otra fila j de la misma matriz
multiplicada por una constante 𝜆.
𝜆𝑅𝑗 + 𝑅𝑖 = 𝑅𝑖
3. Las filas i y j de una matriz pueden ser intercambiadas.
𝑅𝑖 ↔ 𝑅𝑗
El procedimiento para resolver este sistema consta de dos pasos:
1. Eliminación hacia adelante de incógnitas.
2. Sustitución hacia atrás

Ejemplo:
Determinar por el método de eliminación de Gauss el siguiente sistema de ecuaciones:

52
Métodos numéricos

𝑥1 +𝑥2 −𝑥3 = 1 → 𝐸𝑐𝑢𝑎𝑐𝑖ó𝑛 1


3𝑥1 +2𝑥2 +𝑥3 = 1 → 𝐸𝑐𝑢𝑎𝑐𝑖ó𝑛 2
5𝑥1 +3𝑥2 +4𝑥3 = 2 → 𝐸𝑐𝑢𝑎𝑐𝑖ó𝑛 3

Eliminación hacia delante.


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 1 −1 1 3 − 3(1) = 0
[𝟑 𝟐 𝟏 𝟏]Haciendo: Ec.2-3*Ec.1 2 − 3(1) = 2 − 3 = −1
5 3 4 2 1 − 3(−1) = 1+3 =4
{1 − 3(1) = 1 − 3 = −2

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 1 −1 1 5 − 5(1) = 0
[0 −1 4 −2]Haciendo: Ec.3-5*Ec.1 3 − 5(1) = 3 − 5 = −2
𝟓 𝟑 𝟒 𝟐 4 − 5(−1) = 4 + 5 = 9
{2 − 5(1) = 2 − 5 = −3

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 1 −1 1 0 − 2(0) = 0
[0 −1 4 −2]Haciendo: Ec.3-2*Ec.2 −2 − 2(−1) = −2 + 2 = 0
𝟎 −𝟐 𝟗 −𝟑 9 − 2(4) = 9 − 8 = 1
{−3 − 2(−2) = −3 + 4 = 1
Eliminación hacia atrás.

Donde ya se tiene un sistema de ecuaciones en que se pueden despejar las incógnitas.

1 1 −1 1 𝑥1 +𝑥2 −𝑥3 = 1
[𝟎 −1 4 −2] = 0 −𝑥2 +4𝑥3 = −2
𝟎 𝟎 1 1 0 0 𝑥3 = 1
Donde x1=-4, x2=6 y x3=1

Algoritmo para el método de Gauss


function gauss=gauss()
clc;
clear;
d=input('¿Cuantas ecuaciones son?: ');
for i=1:d
for j=1:d
fprintf('Dame el %2.0f coeficiente, de la %2.0f ecuación =
',j,i);
a(i,j)=input(' ');
end
fprintf(' ');
end
for i=1:d
fprintf('Dame el %2.0f termino independiente: ',i);
b(i,1)=input(' ');
end
S=inv(a);
x=S*b;
for i=1:d

53
Métodos numéricos

fprintf(' X%2.0f= %2.2f\n',i, x(i))


end
end

54
Métodos numéricos

Ejercicios
1. Determinar por el método de Eliminación de Gauss los siguientes sistemas de ecuaciones:
1) 𝑥1 +𝑥2 −𝑥3 = −3
6𝑥1 +2𝑥2 +2𝑥3 = 2
−3𝑥1 +4𝑥2 +𝑥3 = 1

2) −𝑥 +5𝑦 +𝑧 = 18
3𝑥 +𝑦 −𝑧 = −2
𝑥 +𝑦 −4𝑧 = −6

3) 3𝑥 −2𝑦 +5𝑧 = 2
3𝑥 −𝑦 −𝑧 = −4
𝑥 +4𝑦 −2𝑧 = −3

4) 5𝑥1 +2𝑥2 +𝑥3 = 3


2𝑥1 +3𝑥2 −3𝑥3 = −10
𝑥1 −3𝑥2 +2𝑥3 = −4

5) 4𝑥1 −𝑥2 +𝑥3 = 8


2𝑥1 +5𝑥2 +2𝑥3 = 3
𝑥1 +2𝑥2 +4𝑥3 = 11

6) 0.2641𝑥1 +0.1735𝑥2 +0.8642𝑥3 = −0.7521


−0.8641𝑥1 −0.4243𝑥2 +0.0711𝑥3 = 0.2501
0.9411𝑥1 +0.0175𝑥2 +0.1463𝑥3 = 0.6310

7) 1 1
𝑥1 + 𝑥2 + 𝑥3 = 2
2 3
1 1 1
𝑥 + 𝑥2 + 𝑥3 = −1
2 1 3 4
1 1 1
𝑥 + 𝑥2 + 𝑥3 = 0
3 1 4 5

2. Una refinería produce gasolina con azufre y sin azufre. Para producir cada tonelada de gasolina
sin azufre se requiere 5 minutos en la planta mezcladora y 4 minutos en la planta de refinación,
mientras que cada tonelada de gasolina con azufre requiere 4 minutos en la planta mezcladora y
2 minutos en la planta de refinación. Si la planta mezcladora está disponible 3 horas y la de
refinación 2 horas, ¿Cuántas toneladas de cada tipo de gasolina deben producirse de modo que
las plantas operen a toda su capacidad? (Kolman y Hill, 2006).

55
Métodos numéricos

3. Tres compuestos se combinan para formar tres tipos de fertilizantes. Una unidad del fertilizante
del tipo I requiere 10 kg del compuesto A, 30 kg del compuesto B y 60 kg del compuesto C. Una
unidad del tipo II requiere 20 kg del A, 30 del B y 50 kg del C. Una unidad del tipo III requiere
50 kg del A y 50 kg del C. Si hay disponibles 1600 kg del A, 1200 kg del B y 3200 del C.
¿Cuántas unidades de los tres tipos de fertilizantes se pueden producir si se usa todo el material
químico disponible?

4. Un fabricante produce reveladores de película de 2, 6 y 9 minutos. La fabricación de cada


tonelada del revelador de 2 minutos requiere 6 minutos en la planta A y 24 minutos en la planta
B. Para manufacturar cada tonelada del revelador de 6 minutos son necesarios 12 minutos en la
planta A y 12 minutos en la planta B. Por último, para producir cada tonelada del revelador de 9
minutos se utiliza 12 minutos la planta A y 12 minutos la planta B. Si la planta A está disponible
10 horas al día y la planta B 16 horas diarias, ¿Cuántas toneladas de cada tipo de revelador de
película pueden producirse de modo que las plantas operen a toda su capacidad? (Kolman y Hill,
2006).

5. Realizar ejercicios: 12.9, 12.18, 12.33 (Chapra y Canale, 2006, p. 336, 342, 345)

56
Métodos numéricos

Método de Gauss-Jordan
El método de Gauss-Jordán es una variante del método de Gauss. Cuando se elimina una
incógnita en una ecuación, Gauss-Jordan elimina esa incógnita en el resto de las
ecuaciones, tomando como base para la eliminación a la ecuación pivote. También todos
los renglones se normalizan cuando se toman como ecuación pivote. El resultado final de
este tipo de eliminación genera una matriz identidad en vez de una triangular como lo
hace Gauss, por lo que no se usa la sustitución hacia atrás para obtener la solución.

Ejemplo:
Resolver el siguiente sistema de ecuaciones con el Método de Gauss-Jordan:

3𝑥1 −0.1𝑥2 −0.2𝑥3 = 7.85


0.1𝑥1 +7𝑥2 −0.3𝑥3 = −19.3
0.3𝑥1 −0.2𝑥2 +10𝑥3 = 71.4

El sistema se expresa como una matriz aumentada

3 −0.1 −0.2 7.85


[0.1 7 −0.3 −19.3]
0.3 −0.2 10 71.4

Ecuación pivote = Ec. 1


Elemento pivote = x1 (incógnita a eliminar de las ecuaciones restantes)

Se normaliza la Ecuación 1

Ec. 1’=Ec.1 (factor) donde factor = (1/3)

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 −0.03333 −0.066667 2.616667 𝐸𝑐. 1′ 3 ∗ (1/3) = 1
[0.1 7 −0.3 −19.3 ] 𝐸𝑐. 2 −0.1 ∗ (1/3) = −0.033333
0.3 −0.2 10 71.4 𝐸𝑐. 3 −0.2 ∗ (1/3) = −0.066667
{ 7.85 ∗ (1/3) = 2.616667

Para obtener la nueva Ec.2:


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
0.1 − (0.1 ∗ 1) = 0
Ec.2= Ec.2 -(0.1)*Ec.1’ 7 − (0.1 ∗ −0.033333) = 7.00333
−0.3 − (0.1 ∗ −0.066667) = −0.29333
{ −19.3 − (0.1 ∗ 2.616667) = −19.5617

Para obtener la nueva Ec.3:

57
Métodos numéricos

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
0.3 − (0.3 ∗ 1) = 0
Ec.3=Ec.3 –(0.3)*Ec.1’ −0.2 − (0.3 ∗ −0.03333) = −0.1900
10 − (0.3 ∗ −0.066667) = −10.0200
{71.4 − (0.3 ∗ 2.616667) = −70.6150

Sistema resultante

1 −0.03333 −0.066667 2.616667 𝐸𝑐. 1


[0 7.00333 −0.29333 −19.5617] 𝐸𝑐. 2
0 −0.19 10.020 70.6150 𝐸𝑐. 3

Ecuación pivote = Ec.2


Elemento pivote = x2(incógnita a eliminar de las ecuaciones restantes)

Se normaliza la Ecuación 2

Ec. 2’=Ec.2 (factor) donde factor = (1/7.00333)


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 −0.03333 −0.066667 2.616667 𝐸𝑐. 1 0 ∗ (1/7.00333) = 0
[0 7.00333 −0.29333 −19.5617] 𝐸𝑐. 2 7.00333 ∗ (1/7.00333) = 1
0 −0.19 10.020 70.6150 𝐸𝑐. 3 −0.29333 ∗ (1/7.00333) = −0.0418848
{ −19.5617 ∗ (1/7.00333) = −2.7932
Se obtiene
1 −0.03333 −0.066667 2.616667 𝐸𝑐. 1
[0 1 −0.0418848 −2.79320] 𝐸𝑐. 2′
0 −0.19 10.020 70.6150 𝐸𝑐. 3

Para obtener la nueva Ec.1:


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 − (−0.03333 ∗ 0) = 1
Ec.1= Ec.1 - (-0.033333)*Ec.2’ −0.03333 − (−0.03333 ∗ 1) = 0
−0.066667 − (0.03333 ∗ −0.0418848) = −0.0680629
{ 2.616667 − (0.03333 ∗ 2.79320) = 2.52356

Para obtener la nueva Ec.3:

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
0 − (−0.19 ∗ 0) = 0
Ec.3=Ec.3 - (-0.19)*Ec.2’ −0.19 − (−0.19 ∗ 1) = 0
10.020 − (−0.19 ∗ −0.0418848) = 10.020
{ 70.6150 − (−0.19 ∗ 2.79320) = 70.0843

Sistema resultante

1 0 −0.0680629 2.52356 𝐸𝑐. 1


[0 1 −0.0418848 −2.79320] 𝐸𝑐. 2
0 0 10.020 70.0843 𝐸𝑐. 3

Ecuación pivote = Ec.3

58
Métodos numéricos

Elemento pivote = x3(incógnita a eliminar de las ecuaciones restantes)

Se normaliza la Ecuación 3

Ec. 3’=Ec.3 (factor) donde factor = (1/10.020)

𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 0 −0.0680629 2.52356 𝐸𝑐. 1 0 ∗ (1/10.020) = 0
[0 1 −0.0418848 −2.79320] 𝐸𝑐. 2 0 ∗ (1/10.020) = 0
0 0 10.020 70.0843 𝐸𝑐. 3 10.020 ∗ (1/10.020) = 1
{70.0843 ∗ (1/10.020) = 7.00003

Se obtiene

1 0 −0.0680629 2.52356 𝐸𝑐. 1


[0 1 −0.0418848 −2.79320] 𝐸𝑐. 2
0 0 1 7.00003 𝐸𝑐. 3′

Para obtener la nueva Ec.1:


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
1 − (−0.0680629 ∗ 0) = 1
Ec.1= Ec.1 - (-0.0680629) *Ec.3’ 0 − (−0.0680629 ∗ 0) = 0
−0.0680629 − (−0.0680629 ∗ 1) = 0
{2.52356 − (−0.0680629 ∗ 7.00003) = 3

Para obtener la nueva Ec.2:


𝑷𝒓𝒐𝒄𝒆𝒅𝒊𝒎𝒊𝒆𝒏𝒕𝒐
0 − (−0.0418848 ∗ 0) = 0
Ec.2= Ec.2 - (-0.0418848) *Ec.3’ 1 − (−0.0418848 ∗ 0) = 1
−0.0418848 − (−0.0418848 ∗ 1) = 0
{−2.79320 − (−0.0418848 ∗ 7.00003) = 2.5

Sistema resultante

1 0 0 3.0 𝐸𝑐. 1
[0 1 0 −2.5] 𝐸𝑐. 2
0 0 1 7.0 𝐸𝑐. 3′

De acuerdo con el resultado, los valores de las incógnitas son:


x1=3, x2=-2.5 y x3=7.0

Algoritmo para el método de Gauss-Jordan


function gauss_jordan=gauss_jordan()
clc;
A= input('Ingrese la matriz: []');
B= input('Ingrese la matriz: []');
C= [A,B];
for i=1: length(C(:,1))
if C(i,i)~=1
C(i,:)=C(i,:)./C(i,i);
disp(C)

59
Métodos numéricos

end
for n=1:length(C(:,1))
if n~=i
C(n,:)=-C(n,i).*C(i,:)+C(n,:);
disp(C)
end
end
end

60
Métodos numéricos

Ejercicios
1. Determinar por el método de Gauss-Jordan los siguientes sistemas de ecuaciones:
1) 10𝑥1 +2𝑥2 −𝑥3 = 27
−3𝑥1 −6𝑥2 +2𝑥3 = −61.5
𝑥1 +𝑥2 +5𝑥3 = −21.5

2) 0.2641𝑥1 +0.1735𝑥2 +0.8642𝑥3 = 0.7521


−0.8641𝑥1 −0.4243𝑥2 −0.0711𝑥3 = 0.2501
0.9411𝑥1 +0.0175𝑥2 +0.1463𝑥3 = 0.6310

3) −0.06𝑥1 +0.04𝑥2 +0.12𝑥3 = 3


0.56𝑥1 −1.56𝑥2 +0.35𝑥3 = 2
−0.24𝑥1 +1.24𝑥2 −0.28𝑥3 = 0

4) 4𝑥1 +𝑥2 +2𝑥4 = 2


𝑥1 −3𝑥2 +𝑥3 = −3
𝑥1 +𝑥2 +4𝑥3 +𝑥4 = 4
𝑥2 +𝑥3 −2𝑥4 = 0

2. Se tienen tres lingotes compuestos del siguiente modo:


• El primero de 20 g de oro, 30 g de plata y 40 g de cobre.
• El segundo de 30 g de oro, 40 g de plata y 50 g de cobre.
• El tercero de 40 g de oro, 50 g de plata y 90 g de cobre.

Se pide, ¿Qué peso habrá de tomarse de cada uno de los lingotes anteriores para formar un nuevo
lingote de 34 g de oro, 46 g de plata y 67 g de cobre?

3. En una acerería se fabrican tres tipos de acero; en láminas, en rollos o aceros especiales. Estos
productos requieren de chatarra, carbón y aleaciones en las cantidades que se indican en la tabla
siguiente por cada unidad de producto fabricado.

Acero en láminas Acero en rollos Acero especial


Chatarra 8 6 6
Carbón 6 6 4
Aleaciones 2 1 3

Si se dispone de 34 unidades de chatarra, 28 de carbón y 9 de aleaciones, ¿Cuántas unidades de


cada tipo de acero se podrán fabricar con estos materiales?

4. Realizar ejercicios: APP3, APP7, APP8 (Gerald y Wheatley ,2004, p.141, p.143, p.144)

61
Métodos numéricos

Métodos Iterativos
Un método iterativo consta de los siguientes pasos.
1. Inicia con una solución aproximada (Semilla).
2. Ejecuta una serie de cálculos para obtener o construir una mejor aproximación
partiendo de la aproximación semilla. La fórmula que permite construir la
aproximación usando otra se conoce como ecuación de recurrencia.
3. Se repite el paso anterior, pero usando como semilla la aproximación obtenida.

Matriz diagonalmente dominante


Una matriz se dice matriz diagonalmente dominante, si en cada uno de los renglones, el
valor absoluto de las incógnitas de la diagonal principal es mayor que la suma de los
valores absolutos de las incógnitas restantes del mismo renglón.
A veces la matriz de un sistema de ecuaciones no es diagonalmente dominante pero cuando
se cambian el orden de las ecuaciones y las incógnitas el nuevo sistema puede tener matriz
de coeficientes diagonalmente dominante.

10 −1 0 4 1 3
[−1 15 −2] 𝑑𝑜𝑚𝑖𝑛𝑎𝑛𝑡𝑒 [2 8 1] 𝑛𝑜 𝑑𝑜𝑚𝑖𝑛𝑎𝑛𝑡𝑒
0 −2 8 3 −10 2

Método de Jacobi o de los desplazamientos simultáneos


Esta técnica muestra cierta similitud con el método de iteración de punto fijo, ya que
consiste en despejar una de las incógnitas de una ecuación dejándola en función de las
otras. La manera más sencilla es despejar a x1 de la primera ecuación, x2 de la segunda
ecuación, xi de la i-ésima ecuación, hasta xn de la n-ésima ecuación

El método Jacobi es el método iterativo para resolver sistemas de ecuaciones lineales más
simple y se aplica solo a sistemas cuadrados, es decir a sistemas con tantas incógnitas como
ecuaciones.

62
Métodos numéricos

Dado el sistema lineal

𝑎1,1 𝑥1 + 𝑎1,2 𝑥2 + … 𝑎1,𝑛 𝑥𝑛 = 𝑏1


𝑎2,1 𝑥1 + 𝑎2,2 𝑥2 + … 𝑎2,𝑛 𝑥𝑛 = 𝑏2
𝑎3,1 𝑥1 + 𝑎3,2 𝑥2 + … 𝑎3,𝑛 𝑥𝑛 = 𝑏3
⋮ ⋮ ⋱ ⋮ ⋮
𝑎𝑚,1 𝑥1 + 𝑎𝑚,2 𝑥2 + … 𝑎𝑚,𝑛 𝑥𝑛 = 𝑏𝑚

Es factible despejar a x1 de la primera ecuación, a x2 de la segunda ecuación, a x3 de la


tercera y así sucesivamente.

𝑥1 = (𝑏1 + 𝑎1,2 𝑥2 +𝑎1,3 𝑥3 + … +𝑎1,𝑛 𝑥𝑛 )/ 𝑎1,1


𝑥2 = (𝑏2 + 𝑎2,1 𝑥1 +𝑎2,3 𝑥3 + … +𝑎2,𝑛 𝑥𝑛 )/ 𝑎2,2
𝑥3 = (𝑏3 + 𝑎3,1 𝑥1 +𝑎3,2 𝑥2 + … +𝑎3,𝑛 𝑥𝑛 )/ 𝑎3,3
⋮ ⋮ ⋮ ⋮ ⋱ ⋮ ⋮
𝑥𝑚 = (𝑏𝑚 + 𝑎𝑚,1 𝑥1 +𝑎𝑚,2 𝑥2 + … +𝑎𝑚,𝑛−1 𝑥𝑛−1 )/ 𝑎𝑚,𝑛

Entonces se parte de una estimación inicial de la solución, x(0), la cual se sustituye en las
ecuaciones para producir una nueva estimación, x(1).
El vector x(1) se sustituye en esas mismas ecuaciones para obtener ahora a x(2). Este
procedimiento se repite entonces para calcular las estimaciones x(3), x(4), x(5), entre otros, y
el proceso termina cuando se cumple alguno de estos criterios de convergencia.
1. ‖𝑥 (𝑘+1) − 𝑥 (𝑘) ‖ < ∈
‖𝑥 (𝑘+1) −𝑥 (𝑘) ‖
2.
‖𝑥 (𝑘+1) ‖
<∈

Uno de los principales problemas de los métodos iterativos es la garantía de que el método
va a converger, es decir, va a producir una sucesión de aproximaciones cada vez
efectivamente más próximas a la solución.

Este método es muy poco utilizado debido a que el método de Gauss-Seidel converge más
rápidamente a la solución y además lo hace cuando no se logra que el método de Jacobi
converja.

La condición suficiente para que el método de Jacobi converja es que la matriz de


coeficientes sea diagonal dominante, es decir que cada elemento de la diagonal principal
es mayor en valor absoluto que la suma del resto de los elementos de la misma fila en la
que se encuentra el elemento en cuestión.

63
Métodos numéricos

Ejemplo:
Con el método de Jacobi, inicie con 𝑥 (0) = 0 y considere 5 decimales como criterio de
convergencia.

Es=(0.5 x 102–n)%
Es=(0.5 x 102-5)% = (0.5 x 10-3)% = 0.0005%= 0.05

10𝑥1 −𝑥2 = 9
−𝑥1 10𝑥2 −2𝑥3 = 7
−2𝑥2 10𝑥3 = 6

Despejando las incógnitas:


x1= (9 + x2 + 0x3 )/10
x2= (7 + x1 + 2x3)/10
x3= (6 +0x1+ 2x2)/10

1ª. Iteración
Con x(0)=[0,0,0] que aplicadas a las estimación inicial x(0) permiten calcular la nueva
iteración x(1).

x1= (9 + 1(0) + 0)/10 = 9/10 =0.9


x2=(7 + 1(0) + 2(0))/10 = 7/10 =0.7
x3=(6 + 0 + 2(0))/10 = 6/10 =0.6
(1) (0)
𝑥 −𝑥 0.9 − 0.0
𝐸𝑎1 = | 1 (1) 1 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥1 0.9
(1) (0)
𝑥 −𝑥 0.7 − 0.0
𝐸𝑎2 = | 2 (1) 2 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥2 0.7

(1) (0)
𝑥 −𝑥 0.6 − 0.0
𝐸𝑎3 = | 3 (1) 3 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥3 0.6

2ª. Iteración
Con x(1)=[0.9, 0.7, 0.6] que aplicadas a las estimación x(1) permiten calcular la nueva
iteración x(2).

x1= (9 + 1(0.7) + 0)/10 = 9.7/10 =0.97


x2=(7 + 1(0.9) + 2(0.6))/10 =9.1 /10 =0.91
x3=(6 + 0 + 2(0.7))/10 = 7.4/10 =0.74
(2) (1)
𝑥 −𝑥 0.97 − 0.9
𝐸𝑎1 = | 1 (2) 1 | 100% = | | 100% = 0.072 ∗ 100% = 7.21
𝑥1 0.97

64
Métodos numéricos

(2) (1)
𝑥 −𝑥 0.91 − 0.7
𝐸𝑎2 = | 2 (2) 2 | 100% = | | 100% = 0.2730 ∗ 100% = 27.30
𝑥2 0.91

(2) (1)
𝑥 −𝑥 0.74 − 0.6
𝐸𝑎3 = | 3 (2) 3 | 100% = | | 100% = 0.1891 ∗ 100% = 18.91
𝑥3 0.74

3ª. Iteración

Con x(2)=[0.97, 0.91, 0.74] que aplicadas a las estimación x(2) permiten calcular la nueva
iteración x(3).

x1= (9 + 1(0.91) + 0)/10 = 9.91/10 =0.9910


x2=(7 + 1(0.97) + 2(0.74))/10 =9.45 /10 =0.9450
x3=(6 + 0 + 2(0.91))/10 = 7.82/10 =0.7820

(3) (2)
𝑥 −𝑥 0.9910 − 0.97
𝐸𝑎1 = | 1 (3) 1 | 100% = | | 100% = 0.0212 ∗ 100% = 2.12
𝑥1 0.9910
(3) (2)
𝑥 −𝑥 0.9450 − 0.91
𝐸𝑎2 = | 2 (3) 2 | 100% = | | 100% = 0.0370 ∗ 100% = 3.70
𝑥2 0.9450

(3) (2)
𝑥 −𝑥 0.7820 − 0.74
𝐸𝑎3 = | 3 (3) 3 | 100% = | | 100% = 0.0537 ∗ 100% = 5.37
𝑥3 0.7820

Nota: Así sucesivamente hasta llegar a 7 iteraciones

iteración x1(k) x2(k) x3(k) Ea1(%) Ea2(%) Ea3(%) Es(%)

semilla 0 0 0 0.05
1 0.9 0.7 0.6 100 100 100
2 0.97 0.91 0.74 7.22 23.08 18.92
3 0.991 0.945 0.782 2.12 3.70 5.37
4 0.9945 0.9555 0.789 0.35 1.10 0.89
5 0.99555 0.95725 0.7911 0.11 0.18 0.27
6 0.99573 0.95778 0.79145 0.018 0.0553 0.0442
7 0.99578 0.95786 0.79156 0.00502 0.0083 0.0138 parar
Nota: Todos los errores son menores al error esperado de 0.05.

Algoritmo para el método de Jacobi


function jacobi=jacobi()
format long;
clc;
e=input('Numero de ecuaciones: ');
A=input('Matriz A [ ]: ');
b=input('Matriz b [ ]: ');
ea=input('Dame el porciento de error: ');
x0=zeros(1,e);

65
Métodos numéricos

err=100;
k=0;
while ea<= err
k=k+1;
fprintf('%2d',k)
for i=1:e
suma=0;
for j=1:e
if i ~=j
suma=suma+A(i,j)*x0(j);
end
end
x(i)=(b(i)-suma)/A(i,i);
fprintf('%10.4f',x(i))
end
err=abs((x-x0)/x)*100;
fprintf('%10.4f\n',err)
x0=x;
end
end

66
Métodos numéricos

Ejercicios
1. Determinar por el método de Jacobi, iniciando con x = 0 y considere 5 decimales como criterio
de convergencia para los siguientes sistemas de ecuaciones.
1) 3.8𝑥 +1.6𝑦 +0.9𝑧 = 15.5
−0.7𝑥 +5.4𝑦 +1.6𝑧 = 10.3
1.5𝑥 +1.1𝑦 −3.2𝑧 = 3.5

2) 31𝑥 +2.2𝑦 +9𝑧 = 82.33


22𝑥 +40𝑦 +2𝑧 = 1112.63
9𝑥 +2𝑦 +31𝑧 = −113.03

3) 6𝑥1 −𝑥2 −𝑥3 +4𝑥4 = 17


𝑥1 −10𝑥2 +2𝑥3 −𝑥4 = −17
3𝑥1 −2𝑥2 +8𝑥3 −𝑥4 = 19
𝑥1 +𝑥2 +𝑥3 −5𝑥4 = −14

4) 5𝑥1 +𝑥2 +2𝑥3 −𝑥4 = 1


𝑥1 +7𝑥2 +3𝑥4 = 2
2𝑥1 +5𝑥3 +𝑥4 = 3
−𝑥1 +3𝑥2 +𝑥3 +8𝑥4 = 4

2. Una compañía de electrónica produce transistores, resistores y chips de computadora.


Cada transistor requiere cuatro unidades de cobre, una de zinc y dos de vidrio. Cada
resistor requiere tres, tres y una unidad de dichos materiales, respectivamente, y cada
chip de computadora requiere dos, una y tres unidades de los materiales,
respectivamente. En forma de tabla, esta información queda así: (Chapra y Canale, 2006,
p. 325)
Componente Cobre Zinc Vidrio
Transistores 4 1 2
Resistores 3 3 1
Chip de computadora 2 1 3

Los suministros de estos materiales varían de una semana a la otra, de modo que la
compañía necesita determinar una corrida de producción diferente cada semana. Por
ejemplo, cierta semana las cantidades disponibles de los materiales son 960 unidades de
cobre, 510 unidades de zinc y 610 unidades de vidrio. Plantee el sistema de ecuaciones
que modela la corrida de producción para resolver cuál es el número de transistores,
resistores y chips de computadora por manufacturar esta semana.

67
Métodos numéricos

3. Determinar las concentraciones molares de una mezcla de cinco componentes en


solución a partir de los siguientes datos espectrofotométricos (Nieves y Domínguez,
2014).

Longitud Absorbancia molar de componente Absorbancia


de onda j total
i 1 2 3 4 5 observada
1 98 9 2 1 0.5 0.1100
2 11 118 9 4 0.88 0.2235
3 27 27 85 8 2 0.2800
4 1 3 17 142 25 0.3000
5 2 4 7 17 118 0.1400

68
Métodos numéricos

Método de Gauss –Seidel


El método de Gauss-Seidel es muy semejante al método de Jacobi. Mientras que en el de
Jacobi se utiliza el valor de las incógnitas para determinar una nueva aproximación, en el de
Gauss-Seidel se va utilizando los valores de las incógnitas recién calculados en la misma
iteración, y no en la siguiente.

Es factible despejar a x1 de la primera ecuación, a x2 de la segunda ecuación, a x3 de la


tercera y así sucesivamente.
(𝑘+1) (𝑘) (𝑘) (𝑘)
𝑥1 = (𝑏1 + 𝑎1,2 𝑥2 +𝑎1,3 𝑥3 + … +𝑎1,𝑛 𝑥𝑛 )/ 𝑎1,1
(𝑘+1) (𝑘) (𝑘) (𝑘)
𝑥2 = (𝑏2 + 𝑎2,1 𝑥1 +𝑎2,3 𝑥3 + … +𝑎2,𝑛 𝑥𝑛 )/ 𝑎2,2
(𝑘+1) (𝑘) (𝑘) (𝑘)
𝑥3 = (𝑏3 + 𝑎3,1 𝑥1 +𝑎3,2 𝑥2 + … +𝑎3,𝑛 𝑥𝑛 )/ 𝑎3,3
⋮ ⋮ ⋮ ⋮ ⋱ ⋮ ⋮
(𝑘+1) (𝑘) (𝑘) (𝑘+1)
𝑥𝑚 = (𝑏𝑚 + 𝑎𝑚,1 𝑥1 +𝑎𝑚,2 𝑥2 + … +𝑎𝑚,𝑛−1 𝑥𝑛−1 )/ 𝑎𝑚,𝑛

Ejemplo:
Con el método de Gauss-Seidel, inicie con x(0)=[0,0,0] y considere 5 decimales como
criterio de convergencia.

Es=(0.5 x 102–n)%
Es=(0.5 x 102-5)% = (0.5 x 10-3)% = 0.0005%= 0.05

10𝑥1 −𝑥2 =9
−𝑥1 10𝑥2 −2𝑥3 =7
−2𝑥2 10𝑥3 =6

Despejando las incógnitas:


x1(k+1)= (9 + 1x2(k) )/10
x2(k+1)= (7 + 1x1(k+1)+ 2x3(k))/10
x3(k+1)= (6 + 2x2(k+1))/10

1ª. Iteración
Con x(0)=[0, 0,0] que aplicadas a las estimación inicial x(0) permiten calcular la nueva
iteración x(1).

x1(1)= (9 + 1x2(0) )/10 = (9 + 1*(0))/10 =0.9


(1) (1) (0)
x2 =(7 + 1x1 + 2x3 )/10 = (7+ 1*(0.9)+ 2*(0))/10 =0.79
x3(1)=(6 + 2x2(1))/10 = (6 + 2*(0.79))/10 =0.758

69
Métodos numéricos

(1) (0)
𝑥 −𝑥 0.9 − 0.0
𝐸𝑎1 = | 1 (1) 1 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥1 0.9
(1) (0)
𝑥 −𝑥 0.79 − 0.0
𝐸𝑎2 = | 2 (1) 2 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥2 0.79

(1) (0)
𝑥 −𝑥 0.758 − 0.0
𝐸𝑎3 = | 3 (1) 3 | 100% = | | 100% = 1 ∗ 100% = 100
𝑥3 0.758

2ª. Iteración
Con x(1)=[0.9, 0.79,0.758] que aplicadas a las estimación inicial x(1) permiten calcular la
nueva iteración x(2).

x1(1)= (9 + 1x2(1) )/10 = (9 + 1*(0.79))/10 =0.979


(1) (2) (1)
x2 =(7 + 1x1 + 2x3 )/10 = (7+ 1*(0.979)+ 2*(0.758))/10 =0.9495
x3(1)=(6 + 2x2(2))/10 = (6 + 2*(0.9495))/10 =0.7899

(2) (1)
𝑥 −𝑥 0.979 − 0.9
𝐸𝑎1 = | 1 (2) 1 | 100% = | | 100% = 0.0806 ∗ 100% = 8.06
𝑥1 0.979
(1) (0)
𝑥 −𝑥 0.9495 − 0.79
𝐸𝑎2 = | 2 (1) 2 | 100% = | | 100% = 0.1679 ∗ 100% = 16.79
𝑥2 0.9495

(1) (0)
𝑥 −𝑥 0.7899 − 0.758
𝐸𝑎3 = | 3 (1) 3 | 100% = | | 100% = 0.0403 ∗ 100% = 4.03
𝑥3 0.7899

3ª. Iteración
Con x(2)=[0.979, 0.9495,0.7899] que aplicadas a las estimación inicial x(2) permiten calcular
la nueva iteración x(3).

x1(2)= (9 + 1x2(2) )/10 = (9 + 1*(0.9495))/10 =0.99495


x2(2)=(7 + 1x1(3) + 2x3(2))/10 = (7+ 1*(0.99495)+ 2*(0.7899))/10 =0.957475
x3(2)=(6 + 2x2(3))/10 = (6 + 2*(0.957475))/10 =0.791495

(2) (1)
𝑥 −𝑥 0.99495 − 0.979
𝐸𝑎1 = | 1 (2) 1 | 100% = | | 100% = 0.01603 ∗ 100% = 1.603
𝑥1 0.99495
(1) (0)
𝑥 −𝑥 0.957475 − 0.9495
𝐸𝑎2 = | 2 (1) 2 | 100% = | | 100% = 0.00832 ∗ 100% = 0.832
𝑥2 0.957475

(1) (0)
𝑥 −𝑥 0.791495 − 0.7899
𝐸𝑎3 = | 3 (1) 3 | 100% = | | 100% = 0.00201 ∗ 100% = 0.201
𝑥3 0.791495

70
Métodos numéricos

Nota: Así sucesivamente hasta llegar a 5 iteraciones

iteración x1(k) x2(k) x3(k) Ea1(%) Ea2(%) Ea3(%) Es(%)


semilla 0 0 0 0.05
1 0.9 0.79 0.758 100 100 100
2 0.979 0.949 0.7899 8.06 16.79 4.03
3 0.995 0.957 0.791 1.603 0.832 0.201
4 0.996 0.958 0.792 0.080 0.042 0.010
5 0.996 0.958 0.792 0.004 0.002 0.001 parar

Algoritmo para el método de Gauss-Seidel


function gaussseidel=gaussseidel()
format long;
clc;
e =input('Numero de ecuaciones: ');
A =input('Matriz A [ ]: ');
b =input('Matriz b [ ]: ');
ea =input('Dame el porciento de error: ');
x0 =zeros(1,e);
x =x0;
err=100;
k=0;
while ea<= err
k=k+1;
fprintf('%2d',k)
for i=1:e
suma=0;
for j=1:e
if i ~=j
suma=suma+A(i,j)*x(j);
end
end
x(i)=(b(i)-suma)/A(i,i);
fprintf('%10.4f',x(i))
end
err=abs((x-x0)/x)*100;
fprintf('%10.4f\n',err)
x0=x;
end
end

71
Métodos numéricos

Ejercicios
1. Determinar por el método de Gauss-Seidel, iniciando con x=0 y considere 5 decimales
como criterio de convergencia para los siguientes sistemas de ecuaciones.

1) 4𝑥1 +0.5𝑥2 +𝑥3 = 8


𝑥1 −10𝑥2 +𝑥3 = −6
−𝑥1 +𝑥2 +5𝑥3 = 10

2) 4𝑥1 +𝑥2 +2𝑥4 = 2


𝑥1 −3𝑥2 +𝑥3 = −3
𝑥1 +𝑥2 +4𝑥3 +𝑥4 = 4
+𝑥2 +𝑥3 −2𝑥4 = 0

3) 10𝑥1 −𝑥2 +2𝑥3 = 6


−𝑥1 +11𝑥2 −𝑥3 +3𝑥4 = 25
2𝑥1 −𝑥2 +10𝑥3 −𝑥4 = −11
+3𝑥2 −𝑥3 +8𝑥4 = 15

4) 4.63𝑥1 −1.21𝑥2 +3.22𝑥3 = 2.22


−3.07𝑥1 +5.48𝑥2 +2.11𝑥3 = −3.17
1. 26𝑥1 +3.11𝑥2 +4.57𝑥3 = 5.11

2. Un ingeniero civil que trabaja en la construcción requiere 4800, 5800 y 5700 m3 de arena,
grava fina, y grava gruesa, respectivamente, para cierto proyecto constructivo. Hay tres
canteras de las que puede obtener dichos materiales. La composición de dichas canteras es
la que sigue. (Chapra y Canale, 2006, p. 342)

Arena % Grava fina % Grava gruesa %


Cantera 1 55 30 15
Cantera 2 25 45 30
Cantera 3 25 20 55

¿Cuántos metros cúbicos se deben extraerse de cada cantera a fin de satisfacer las
necesidades del ingeniero?

3. En una acerería se fabrican tres tipos de acero; en láminas, en rollos o aceros especiales.
Estos productos requieren de chatarra, carbón y aleaciones en las cantidades que se indican
en la tabla siguiente por cada unidad de producto fabricado.

Acero en láminas Acero en rollos Acero especial


Chatarra 8 2 5
Carbón 5 6 1
Aleaciones 2 1 7

Si para el próximo mes el ingeniero encargado de su fabricación requiere 170 unidades de


acero en láminas, 100 unidades de acero en rollos y 110 unidades de acero especial, indique
las unidades de chatarra, carbón y aleaciones que se requieren.

4. Un albañil compra cemento, cal y arena en la ferretería más cercana

72
Métodos numéricos

En la primera vez, compro 6 bultos de cemento, 2 bultos de cal y un viaje de arena y pago
$2,064 pesos.
En la segunda vez compro 2 bultos de cemento, 4 de cal y medio viaje de arena, pagando
$1,028 pesos, pero no le alcanzo el cemento ni la arena, así que compro 1 bulto de cemento
y 2 viajes de arena pagando $2,160 pesos.
El albañil quiere saber cual es el precio por bulto de cemento, cal y por el viaje de arena.

5. Para la producción de cuatro tipos de computadoras, se requieren cuatro clases de recursos,


horas/hombre, metales, plásticos y componentes electrónicos, en la producción. En el
cuadro se resume las cantidades necesarias para cada uno de estos recursos en la producción
de cada tipo de computadora. Si se dispone diariamente de 504 horas/hombre, 412 kg de
metal, 297 kg de plástico y 350 componentes electrónicos, ¿Cuántas computadoras de cada
tipo se pueden construir por día?

Horas/hombre Metales Plásticos Componentes


Computadora
Kg/computadora Kg/computadora Kg/computadora Unidades/computadora
1 23 9 6 4
2 4 21 5 8
3 7 6 35 5
4 12 7 9 32

5. Realizar ejercicios: 3.63 (Nieves y Domínguez, 2006, p. 253)

73
Métodos numéricos

Tabla 1. Ventajas y desventajas de los métodos iterativos comparados con los métodos directos
(Nieves y Domínguez, 2014, p. 247)
Ventajas Desventajas
1. Probablemente más eficientes que los 1. Si se tienen varios sistemas que comparten la
directos para sistemas de orden muy alto. matriz coeficiente, esto no representará ahorro
2. Más simples de programar. de cálculos ni de tiempo de máquina, ya que
3. Puede aprovecharse una aproximación a la por cada vector a la derecha de A tendrá que
solución, si tal aproximación existe. aplicarse el método seleccionado.
4. Se obtienen aproximaciones burdas de la 2. Aun cuando la convergencia esté asegurada,
solución con facilidad. puede ser lenta y, por tanto, los cálculos
5. Son menos sensibles a los errores de requeridos para obtener una solución
redondeo (valioso en sistemas mal particular no son predecibles.
condicionados). 3. El tiempo de máquina y la exactitud del
6. Se requiere menos memoria de máquina. resultado dependen del criterio de
Generalmente, las necesidades de memoria convergencia.
son proporcionales al orden de la matriz. 4. Si la convergencia es lenta, los resultados
deben interpretarse con cautela.
5. No se tiene ventaja particular alguna (tiempo
de máquina por iteración) si la matriz
coeficiente es simétrica.

74
Métodos numéricos

INTERPOLACIÓN
Interpolación
La interpolación es, a partir de una seria de puntos, obtener una ecuación cuya curva pase
por todos ellos o lo más cerca posible.
La extrapolación es el proceso de calcular un valor de 𝑓(𝑥) que cae fuera del rango de los
puntos base conocidos x0, x1,... , xn. La interpolación más exacta usualmente se obtiene
cuando las incógnitas caen cerca de los puntos base.
Obviamente, esto no sucede cuando las incógnitas caen fuera del rango y, por lo tanto, el
error en la extrapolación puede ser muy grande. La naturaleza abierta en los extremos de la
extrapolación representa un paso en la incógnita porque el proceso extiende la curva más
allá de la región conocida. Como tal, la curva verdadera diverge fácilmente de la
predicción. Por lo tanto, se debe tener cuidado extremo en casos donde se deba extrapolar.

La idea básica de la interpolación es hallar un polinomio o función que cumpla con pasar
por todos los puntos de un conjunto de datos (x1, y1), (x2, y2),…,(xn, yn), y poder estimar los
valores entre ellos por medio del polinomio.
El método más común para interpolar valores intermedios es la interpolación polinomial, la
cual consiste en determinar el polinomio de orden n que ajusta a n+1 datos. La
interpolación de Lagrange es una de las alternativas más atractivas que existe para
interpolar, debido a la facilidad de programar.

Para generar una interpolación de orden n, es necesario contar con n+1 datos conocidos,
por ejemplo:
• Una interpolación lineal requiere de conocer dos datos.
• Una interpolación cuadrática requiere de conocer tres datos.
• Una interpolación cúbica requiere de conocer cuatro datos.

75
Métodos numéricos

Aproximación polinomial simple e interpolación


La interpolación es de gran importancia en el campo de la ingeniería, ya que, al consultar
fuentes de información presentadas en forma tabular, es frecuente no encontrar el valor
buscado como un punto en la tabla, por ejemplo, las tablas 1 y 2 presentan la temperatura
de ebullición de la acetona (𝐶3 𝐻6 𝑂) a diferentes presiones.
Tabla 1
Puntos 0 1 2 3 4 5 6
P (atm) 1 2 5 10 20 30 40
0
T ( C) 56.5 78.6 113.0 144.5 181.0 205.0 214.5

Tabla 2
Puntos 0 1 2 3
P (atm) 1 5 20 40
T (0C) 56.5 113.0 181.0 214.5

Supóngase que solo se dispusiera de la tabla 2 y se desea calcular la temperatura de


ebullición de la acetona a 2 atm de presión.
Una forma muy común de resolver este problema es sustituir los puntos (0) y (1) de la tabla
en la ecuación de la línea recta: P(x)=a0+a1x, de tal modo que resultan 2 ecuaciones con 2
incógnitas que son a0 y a1. Con la solución del sistema se consigue una aproximación
polinomial de primer grado, lo que permite efectuar interpolaciones lineales; es decir:

Sistema de ecuaciones para 2 puntos (0) y (1).


56.5 =a0 + a1x
113 =a0 + a1x

Sustituyendo los valores de P(atm) del punto (0) y punto (1) en las x del sistema de
ecuaciones para obtener el valor de T(0C) se tiene:

56.5 = a0 + a1(1) para el punto (0)


113 = a0 + a1(5) para el punto (1)

Sistema que al resolver (por sustitución) nos da: a0=42.375 y a1=14.125


Por tanto, estos valores generan la ecuación: P1(x)=42.375 + 14.125x

76
Métodos numéricos

La ecuación resultante puede emplearse para aproximar la temperatura cuando la presión es


conocida.

P(2)=42.375 +14.125(2)=70.625 0C

Al sustituir la presión x=2 atm se obtiene la temperatura de 70.6 0C.

Figura 6. Interpolación de la acetona con 2 nodos (Nieves y Domínguez, 2006)

A este proceso se le conoce como interpolación, si se quisiera una aproximación mejor al


valor “verdadero” de la temperatura buscada, podrían unirse más puntos de la tabla con una
curva suave (sin picos), por ejemplo, tres puntos (0), (1), (2).
Analíticamente el problema se resuelve al aproximar la función desconocida T=f(P) con un
polinomio que pasa por los tres puntos (0), (1), (2).

Continuando con el ejemplo anterior para tres puntos sustituyendo en las x los valores:

56.5=a0 + a1(1) + a2(1)2


113=a0 + a1(5) + a2(5)2
181=a0 + a1(20) + a2(20)2

Al resolver el sistema de ecuaciones (resolver con Gauss-Jordan) se obtienen los siguientes


resultados:

77
Métodos numéricos

a0=39.8509 a1=17.1539 a2=-0.5048


De tal modo que la ecuación polinomial queda: P2(x)=39.8509 +17.1539x -0.5048x2
Calculando nuevamente para x=2 atm en la nueva ecuación nos da:

P(2)= 39.8509 +17.1539(2) -0.5048(2)2=72.1395 0C.

Figura 7. Interpolación de la acetona con 3 nodos (Nieves y Domínguez, 2006)

La aproximación a la temperatura “correcta” es obviamente mejor en este caso. Obsérvese


que ahora se ha aproximado la función desconocida T=f(P) con un polinomio de segundo
grado (parábola) que pasa por los tres puntos más cercanos al valor buscado. En general, si
se desea aproximar una función con un polinomio de grado n, se necesitan más puntos que
sustituidos en la ecuación polinomial de grado n.

Continuar el ejemplo, agregando un nuevo nodo para resolver el sistema de ecuaciones,


indicar cuál es su polinomio, su gráfica, así como el resultado.

Una vez resuelto el sistema se sustituyen los valores de ai en la ecuación polinomial con el
cual se obtiene el polinomio de aproximación. A este método se le conoce como
aproximación polinomial simple y tiene la forma general:

𝑷𝒏 (𝒙) = 𝒂𝟎 + 𝒂𝟏 𝒙 + 𝒂𝟐 𝒙𝟐 + ⋯ + 𝒂𝒏 𝒙𝒏 [𝑃𝑜𝑙𝑖𝑛𝑜𝑚𝑖𝑎𝑙 𝑑𝑒 𝑔𝑟𝑎𝑑𝑜 𝑛]

También se pueden tener funciones conocidas pero complicadas, por ejemplo:

78
Métodos numéricos

2 1/2
𝑓(𝑥) = ( ) 𝑆𝑒𝑛 𝑥
𝑥
Lo cual conviene, para propósitos prácticos, aproximar con otra función más sencilla, como
un polinomio. El procedimiento es generar una tabla de valores mediante la función
original y a partir de dicha tabla aplicar el método descrito anteriormente.
Algoritmo de la aproximación polinomial simple
function aproximacionpolinomialsimple=aproximacionpolinomialsimple()
format long;
clear X;
clear Y;
clc;
X=input('Dame el vector de puntos en X [ ]: ');
Y=input('Dame el vector de puntos en Y [ ]: ');
x=input('Dame el punto a interpolar: ');
xx=min(X):0.01:max(X);
yy=spline(X,Y,xx);
hold on
plot(X,Y,'o',xx,yy);
grid,title('\bf Aproximación Polinomial Simple')
y=0;
w=length(X);
for i=1:w-1
if x>=X(i) && x<=X(i+1)
y=(Y(i+1)-Y(i))/(X(i+1)-X(i))*(x-X(i))+Y(i);
end
end
fprintf('\n\nResultado de la interpolación lineal es=%12.6f \n',y);
plot(x,y,'sr');
hold off
end

Solución en Matlab:
Dame el vector de puntos en X [ ]: [1, 5, 20]
Dame el vector de puntos en Y [ ]: [56.5, 113, 181]
Dame el punto a interpolar: 2

Resultado de la interpolación lineal es= 70.625000

79
Métodos numéricos

Ejercicios
1) La población de la ciudad de Xalapa está dada por la siguiente tabla:

Nodo 0 1 2 3 4 5
año 1940 1950 1960 1970 1980 1990
habitantes 132,165 151,326 179,323 203,302 226,542 249,633

a) Obtener un polinomio de grado superior a 2 que, según su juicio, sea la mejor aproximación
de los datos originales. Mostrar gráfica.
b) Obtener la población en los años 1948, 1951, 1968, 1975 y 1988.

2) Se tiene la función f(x)= ln (x+1)


a) Evaluar los puntos f(0), f(0.6) y f(0.9)
b) Encontrar un polinomio pn(x) de grado 2 que aproxime al punto f(0.45)
c) Mostrar gráficamente los resultados.

3) A continuación, se presentan las presiones de vapor del cloruro de magnesio.


Nodo 0 1 2 3 4 5 6 7
P(mmHg) 10 20 40 63 100 200 400 760
T(°C) 930 988 1050 1088 1142 1316 1223 1418
Calcule la presión de vapor correspondiente a T=1000 °C.

4) En estudios de polimerización inducida por radiación, se empleó una fuente de rayos gamma
para administrar dosis medidas de radiación. Sin embargo, la dosificación varió con la posición
en el aparato, con estas cifras.

Nodo 0 1 2 3 4 5 6 7
Posición 0 0.5 1.0 1.5 2.0 3.0 3.5 4.0
5
Dosis 10 , rads/hrs. 1.90 2.39 2.71 2.98 3.20 3.20 2.98 2.74

Por alguna razón, no se informó la lectura a 2.5 pulgadas, pero se necesita el valor de la radiación
allí. Ajuste los polinomios interpoladores de varios grados a los datos para proporcionar la
información que falta. Qué hacer ¿Crees que es la mejor estimación para el nivel de dosificación de
2.5 pulgadas?

80
Métodos numéricos

Polinomios de Lagrange
El método de aproximación polinomial simple requiere la solución de un sistema de
ecuaciones algebraicas lineales que, cuando el grado del polinomio es alto, puede presentar
inconvenientes. Existen otros métodos de aproximación polinomial en que no se requiere
resolver un sistema de ecuaciones lineales y los cálculos se realizan directamente; entre
estos se encuentra el de aproximación polinomial de Lagrange.

Se parte nuevamente de una función desconocida f(x) dada en forma tabular y se asume
que:

Para un polinomio de primer grado (ecuación de la recta) puede escribirse:


(𝑥 − 𝑥1 ) (𝑥 − 𝑥0 )
𝑃1 (𝑥) = 𝑦0 + 𝑦1
(𝑥0 − 𝑥1 ) (𝑥1 − 𝑥0 )

Para un polinomio de segundo grado:


(𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥2 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥1 )
𝑃2 (𝑥) = 𝑦0 + 𝑦1 + 𝑦2
(𝑥0 − 𝑥1 )(𝑥0 − 𝑥2 ) (𝑥1 − 𝑥0 )(𝑥1 − 𝑥2 ) (𝑥2 − 𝑥0 )(𝑥2 − 𝑥1 )

Para un polinomio de tercer grado:


(𝑥 − 𝑥1 )(𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥2 )(𝑥 − 𝑥3 )
𝑃3 (𝑥) = 𝑦0 + 𝑦1
(𝑥0 − 𝑥1 )(𝑥0 − 𝑥2 )(𝑥0 − 𝑥3 ) (𝑥1 − 𝑥0 )(𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 )
(𝑥 − 𝑥0 )(𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥1 )(𝑥 − 𝑥2 )
+ 𝑦2 + 𝑦3
(𝑥2 − 𝑥0 )(𝑥2 − 𝑥1 )(𝑥2 − 𝑥3 ) (𝑥3 − 𝑥0 )(𝑥3 − 𝑥1 )(𝑥3 − 𝑥2 )

Por inducción para obtener polinomios de cuarto o n-ésimo grado; este último queda como
se indica a continuación.

(𝑥 − 𝑥0 ) … (𝑥 − 𝑥𝑘−1 )(𝑥 − 𝑥𝑘+1 ) … (𝑥 − 𝑥𝑛 )


𝐿𝑛,𝑘 (𝑥) =
(𝑥𝑘 − 𝑥0 ) … (𝑥𝑘 − 𝑥𝑘−1 )(𝑥𝑘 − 𝑥𝑘+1 ) … (𝑥𝑘 − 𝑥𝑛 )

Que en forma más compacta queda de la siguiente forma:


𝒏

𝑷𝒏 (𝒙) = ∑ 𝒚𝒌 𝑳𝒏,𝒌 (𝒙) [𝑃𝑜𝑙𝑖𝑛𝑜𝑚𝑖𝑜 𝑑𝑒 𝐿𝑎𝑔𝑟𝑎𝑛𝑔𝑒]


𝒌=𝟎

Ejemplo:
a) Obtener la aproximación polinomial de Lagrange con todos los puntos.
b) Interpole el valor de la función f(x) para x=1.8.

81
Métodos numéricos

puntos 0 1 2 3
xi 0 1 3 6
f(x) -3 0 5 7

a) Observe que se tienen 4 puntos en la tabla, por lo que el polinomio será de tercer
grado que se sustituyen en la ecuación.

Sustituimos los valores

(𝑥 − 𝑥1 )(𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥2 )(𝑥 − 𝑥3 )


𝑃3 (𝑥) = 𝑦0 + 𝑦1
(𝑥0 − 𝑥1 )(𝑥0 − 𝑥2 )(𝑥0 − 𝑥3 ) (𝑥1 − 𝑥0 )(𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 )
(𝑥 − 𝑥0 )(𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥0 )(𝑥 − 𝑥1 )(𝑥 − 𝑥2 )
+ 𝑦2 + 𝑦3
(𝑥2 − 𝑥0 )(𝑥2 − 𝑥1 )(𝑥2 − 𝑥3 ) (𝑥3 − 𝑥0 )(𝑥3 − 𝑥1 )(𝑥3 − 𝑥2 )

(𝑥 − 1)(𝑥 − 3)(𝑥 − 6) (𝑥 − 0)(𝑥 − 3)(𝑥 − 6) (𝑥 − 0)(𝑥 − 1)(𝑥 − 6)


𝑃3 (𝑥) = (−3) + (0) + (5)
(0 − 1)(0 − 3)(0 − 6) (1 − 0)(1 − 3)(1 − 6) (3 − 0)(3 − 1)(3 − 6)
(𝑥 − 0)(𝑥 − 1)(𝑥 − 3)
+ (7)
(6 − 0)(6 − 1)(6 − 3)

(𝑥 − 1)(𝑥 − 3)(𝑥 − 6) (𝑥 − 0)(𝑥 − 3)(𝑥 − 6) (𝑥 − 0)(𝑥 − 1)(𝑥 − 6)


𝑃3 (𝑥) = (−3) + (0) + (5)
−18 10 −18
(𝑥 − 0)(𝑥 − 1)(𝑥 − 3)
+ (7)
90

Efectuamos las operaciones

1 (𝑥 3 − 10𝑥 2 + 27𝑥 − 18) −5 (𝑥 3 − 7𝑥 2 + 6𝑥) 7 (𝑥 3 − 4𝑥 2 + 3𝑥)


𝑃3 (𝑥) = ( ) + ( ) +( )
6 18 90

Y finalmente resulta:

1 3 1 2 276
𝑃3 (𝑥) = − 𝑥 − 𝑥 + 𝑥−3
30 30 90

b) El valor de x=1.8 se sustituye en la aproximación polinomial de Lagranje de tercer


grado obtenida anteriormente y el resultado es el siguiente.

1 1 276
𝑃3 (1.8) = − (1.8)3 − (1.8)2 + (1.8) − 3 = 𝟐. 𝟐𝟏
30 30 90

82
Métodos numéricos

Algoritmo para el polinomio de Lagrange


function intlagrange=intlagrange()
format short;
clear x;
clear y;
clc;
X=input('Dame el vector de puntos en X [ ]: ');
Y=input('Dame el vector de puntos en Y [ ]: ');
x =input('Dame el punto a interpolar: ');
n1=length(X);
n=n1-1;
M=zeros(n1,n1);
xx=min(X):0.01:max(X);
yy=spline(X,Y,xx);
hold on
plot(X,Y,'o',xx,yy);
grid,title('\bf Interpolación con el Método de Lagrange')
y=0;
w=length(X);
for i=1:w
L=1;
for j=1:w
if j~=i
L=L*(x-X(j))/(X(i)-X(j));
end
end
y=y+L*Y(i);
end
for k=1:n+1
v=1;
for j=1:n+1;
if k~=j
v=conv(v,poly(X(j)))/(X(k)-X(j));
end
end
M(k,:)=v;
C=Y*M;
Pol=C;
end
fprintf('\n\nResultado de la interpolación es=%12.6f \n',y);
fprintf('\n\nEl Polinomio es de grado=%12.6f \n',n);
disp(Pol)
plot(x,y,'sr');
hold off
end

Solución en Matlab:

83
Métodos numéricos

84
Métodos numéricos

Ejercicios
1. Realice el siguiente ejercicio.
a. Obtener la aproximación polinomial de Lagrange con todos los puntos.
b. Interpole el valor del logaritmo natural de 3.

puntos 0 1 2 3
x 1 4 6 8
ln x 0 1.386294 1.791759 2.079441

2. Determine la ecuación que vincula a los siguientes datos.

puntos 0 1 2 3 4 5 6
x 129 247 530 1550 3010 4820 8010
y 9.46 8.28 5.26 2.77 2.16 1.98 1.22

3. Encuentre los valores de la variable dependiente para x=4.125, 4.375, 5.896, 9.788, 10.500,
10.788, 10.987, dado los siguientes datos.

puntos 0 1 2 3 4 5 6 7
x 2.156 3.145 6.725 7.222 8.434 9.525 10.112 11.028
y 8.112 12.322 14.580 12.366 24.845 28.366 30.554 38.687

4. En la siguiente tabla se presentan los alargamientos de un resorte correspondiente a fuerzas


de diferentes magnitudes que lo deforman.
Puntos 0 1 2 3 4
Fuerza (kgf): x 0 2 3 6 7
Longitud del resorte (m): y 0.120 0.153 0.170 0.225 0.260

Obtener la Interpolación de p (1).

5. Al medir la velocidad (con un tubo Pitot) en una tubería circular de diámetro interior de 20
cm, se encontró la siguiente información:
v (cm/s) 600 550 450 312 240
r (cm) 0 3 5 7 8
Donde r es la distancia en cm medida a partir del centro del tubo, Calcule la velocidad en el
punto r= 4 cm.

85
Métodos numéricos

Polinomio de Newton (Diferencias divididas)


La forma de Lagrange es útil si con los mismos nodos se prueban otros valores de la
función u otras funciones pues las modificaciones a realizar no afectan a los polinomios de
lagrange. En cambio, si se aumenta un nodo los cálculos hay que rehacerlos desde el
principio. Para esta situación es más conveniente la forma de Newton del polinomio de
interpolación.

Uno de los inconvenientes de los polinomios interpoladores de Lagrange es que no hay


relación entre la construcción de Pn-1(x) y la de Pn(x); cada polinomio debe construirse
individualmente y se requieren muchas operaciones para calcular polinomios de grado
elevado.
Aplicando la recurrencia a los polinomios pn, pn-1,…,p1, p0 con p0(x) = f(x0) = f[x0], el
polinomio de interpolación se puede expresar de la llamada forma de Newton como:

𝑷(𝒙)
= 𝒇[𝒙𝟎 ] + 𝒇[𝒙𝟎 , 𝒙𝟏 ](𝒙 − 𝒙𝟎 ) + 𝒇[𝒙𝟎 , 𝒙𝟏 , 𝒙𝟐 ](𝒙 − 𝒙𝟎 )(𝒙 − 𝒙𝟏 ) + ⋯
+ (𝒙 − 𝒙𝟎 )(𝒙 − 𝒙𝟏 ) … 𝒇[𝒙𝟎 , 𝒙𝟏 , … 𝒙𝒏 )(𝒙
− 𝒙𝒏−𝟏 ) [𝑃𝑜𝑙𝑖𝑛𝑜𝑚𝑖𝑜 𝑑𝑒 𝑁𝑒𝑤𝑡𝑜𝑛]

Tabla para calcular el polinomio de Newton (diferencias divididas)


Orden
x f(x)
Primer (a) Segundo (b) Tercer (c) Cuarto (d) n
𝑥0 𝑓[𝑥0 ]
𝑓(𝑥1 ) − 𝑓(𝑥0 )
𝑎0 =
(𝑥1 − 𝑥0 )
𝑎1 − 𝑎0
𝑥1 𝑓[𝑥1 ] 𝑏0 =
(𝑥2 − 𝑥0 )
𝑓(𝑥2 ) − 𝑓(𝑥1 ) 𝑏1 − 𝑏0
𝑎1 = 𝑐0 =
(𝑥2 − 𝑥1 ) (𝑥3 − 𝑥0 )
𝑎2 − 𝑎1 𝑐1 − 𝑐0
𝑥2 𝑓[𝑥2 ] 𝑏1 = 𝑑0 =
(𝑥3 − 𝑥1 ) (𝑥4 − 𝑥0 )
𝑓(𝑥3 ) − 𝑓(𝑥2 ) 𝑏2 − 𝑏1
𝑎2 = 𝑐1 =
(𝑥3 − 𝑥2 ) (𝑥4 − 𝑥1 )
𝑎3 − 𝑎2
𝑥3 𝑓[𝑥3 ] 𝑏2 =
(𝑥4 − 𝑥2 )
𝑓(𝑥4 ) − 𝑓(𝑥3 )
𝑎3 =
(𝑥4 − 𝑥3 )
𝑥4 𝑓[𝑥4 ]

Ejemplo:
La información de la siguiente tabla se obtuvo de un polinomio

86
Métodos numéricos

Puntos 0 1 2 3 4 5
x -2 -1 0 2 3 6
f(x) -18 -5 -2 -2 7 142
A partir de ella:

a) Elabore una tabla de diferencias divididas que interpole la tabla dada.


b) Obtener la interpolación de p (1).

Solución:

a) Tabla del polinomio de Newton.

Polinomio de Newton (Diferencias divididas)


x f(x) Orden
Primer (a) Segundo (b) Tercer (c) Cuarto (d) Quinto (e)
𝑥0 -2 𝑓[𝑥0 ] -18
−5 − (−18)
𝑎0 = = 13
−1 − (−2)
3 − 13
𝑥1 -1 𝑓[𝑥1 ] -5 𝑏0 = = −5
0 − (−2)
−2 − (−5) −1 − (−5)
𝑎1 = =3 𝑐0 = =1
0 − (−1) 2 − (−2)
0−3 1−1
𝑥2 0 𝑓[𝑥2 ] -2 𝑏1 = = −1 𝑑0 = =0
2 − (−1) 3 − (−2)
−2 − (−2) 3 − (−1) 0−0
𝑎2 = =0 𝑐1 = =1 𝑒0 = =0
2−0 3 − (−1) 6 − (−2)
9−0 1−1
𝑥3 2 𝑓[𝑥3 ] -2 𝑏2 = =3 𝑑1 = =0
3−0 6 − (−1)
7 − (−2) 9−3
𝑎3 = =9 𝑐2 = =1
3−2 6−0
45 − 9
𝑥4 3 𝑓[𝑥4 ] 7 𝑏3 = =9
6−2
142 − 7
𝑎4 = = 45
6−3
𝑥5 6 𝑓[𝑥5 ] 142

Polinomio de Newton (Diferencias divididas)


Orden
punto x f(x)
1er. 2o. 3er. 4º. 5º.
0 -2 -18
13
1 -1 -5 -5
3 1
2 0 -2 -1 0
0 1 0
3 2 -2 3 0
9 1
4 3 7 9
45
5 6 142

87
Métodos numéricos

Dónde:
𝑃5 (𝑥) = 𝑓[𝑥0 ] + 𝑓[𝑥0 , 𝑥1 ](𝑥 − 𝑥0 ) + 𝑓[𝑥0 , 𝑥1 , 𝑥2 ](𝑥 − 𝑥0 )(𝑥 − 𝑥1 ) + ⋯
+ (𝑥 − 𝑥0 )(𝑥 − 𝑥1 ) … 𝑓[𝑥0 , 𝑥1 , … 𝑥𝑛 )(𝑥 − 𝑥𝑛−1 )

P5(x)=-18 +13(x-(-2)) -5(x-(-2))(x-(-1)) +1(x-(-2))(x-(-1))(x-0) + 0(x-(-2))(x-(-1))(x-0)(x-2)


+ 0(x-(-2))(x-(-1))(x-0)(x-2)(x-3) =

P5(x)=-18 +13(x +2) -5(x +2)(x + 1) +(x +2)(x +1)(x-0) =

P5(x)=-18 +13x +26 -5(x2 + 3x + 2) +(x3 + 3x2 +2x) =

P5(x)=-18 +13x +26 -5x2 - 15x - 10 + x3 + 3x2 +2x =

P5(x)= -2 - 2x2 + x3

b) Interpolación de p(1).
P5(1)=-2 -2(1)2+ (1)3 = -2-2+1= -4 +1 = -3

Algoritmo pare el polinomio de Newton (diferencias divididas)


function newton=newton()
clc;
clear;
x=input('Dame el vector de puntos en X [ ]: ');
y=input('Dame el vector de puntos en Y [ ]: ');
xu=x;
yu=y;
%Formamos una Matriz cero de la longitud del Vector "Y"
d=zeros(length(y));
% Se genera una Matriz cero donde la primera columna sea el Vector "Y"
d(:,1)=y';
%Formación de las diferencias divididas
for k=2:length(x)
for j=1:length(x)+1-k
d(j,k)=(d(j+1,k-1)-d(j,k-1))/(x(j+k-1)-x(j));
end
end
% Formación del Polinomio
for w=1:length(x)
dq=num2str(abs(d(1,w)));
if w>1
if x(w-1)<0
signo1=' +';
else
signo1=' -';
end
end
if d(1,w)<0
signo2=' -';
else
signo2=' +';
end
if w==1

88
Métodos numéricos

acum=num2str(d(1,1));
elseif w==2
polinomioactual=['(x' signo1 num2str(abs(x(w-1))) ')' ];
actual=[dq '*' polinomioactual];
acum=[acum signo2 actual];
else
polinomioactual=[polinomioactual '.*' '(x' signo1
num2str(abs(x(w-1))) ')' ];
actual=[dq '*' polinomioactual];
acum=[acum signo2 actual];
end
end
% Presentación de Resultados
fprintf('\n Valores de X y Y \n');
disp(xu);
disp(yu);
fprintf('\n Polinomio de Interpolación : %s \n\n',acum);
x=input(' X interp = ');
if x>max(xu)||x<min(xu)
fprintf('\t Valor de "X" fuera de rango. El resultado puede ser
equivocado \n\a');
end
xinterp=x;
yinterp=eval(acum);
fprintf(' Y interp (%g) = %g \n',x,yinterp);
% Grafica de los puntos
fprintf('\n Presiona una tecla para graficar los datos\n');
pause
xg=linspace(min(xu),max(xu));
x=xg;
yg=eval(acum);
plot(xg,yg,xu,yu,'xk',xinterp,yinterp,'sr');
grid,title('\bf Diferencias Divididas de Newton')
end

89
Métodos numéricos

Ejercicios
1. Aplicando interpolación del polinomio de Newton encuentre, para los datos que se
presentan en la tabla que sigue, los valores de la variable dependiente para x= 4.3, 4.5, 5.8 y
9.7
puntos 0 1 2 3 4 5 6
x 2.2 4.2 6.2 8.2 10.2 12.2 14.2
y 10 15 14 16 17 19 20

2. Encuentre el valor de y para x=22.4 y 32.5 dados los siguientes datos.

puntos 0 1 2 3 4 5 6
x 22.2 24.2 26.5 28.3 30.2 32.7 34.0
y 100 150 140 160 170 190 200

3. Obtener una aproximación de f(2.1) usando todos los datos de la siguiente tabla.

puntos 0 1 2 3 4
x 2.0 2.2 2.4 2.6 2.8
f(x) 0.5103757 0.5207843 0.5104147 0.4813306 0.4359160

4. Una de las propiedades que comúnmente se emplean en mecánica de fluidos es el volumen


del líquido. El volumen de un líquido es una función de la temperatura. El líquido más
utilizado por el hombre es el H2O. A continuación, se muestra el volumen de un gramo de
H2O, en el intervalo de 273.150 K a 279.150 K.

puntos 0 1 2 3 4 5 6
0
Temperatura ( K) 273.15 274.15 275.15 276.16 277.15 278.15 279.15
1.000132
Volumen 1.0000733 1.0000321 1.0000078 1.0000000 1.0000081 1.0000318
9

Determinar el volumen para las temperaturas de: 274, 275, 277, 278, 279.

5. Considere la siguiente tabla de salarios:

Salario (Pesos) 0-1,000 1,000-2,000 2,000-3,000 3,000-4,000


Frecuencia 9 30 35 42
Estimar la cantidad de personas con salarios entre 1,000 y 1,500 pesos

90
Métodos numéricos

INTEGRACIÓN NUMÉRICA
De acuerdo con la definición del diccionario, integrar significa "unir todas las partes en un
todo"; unificar; indicar la cantidad total, suma total...".
Matemáticamente, la integración se representa por:

𝑏
𝐼 = ∫ 𝑓(𝑥)𝑑𝑥
𝑎

La cual representa a la integral de la función f(x) con respecto a la variable x, evaluada


entre los límites x=a y x=b.
Como lo sugiere la definición del diccionario, el "significado" de la ecuación es el valor
total o sumatoria de f(x) dx sobre el intervalo desde x=a hasta b.
En realidad, el símbolo∫ es una S mayúscula estilizada que indica la conexión cercana
entre la integración y la sumatoria.
Para las funciones que se encuentran sobre el eje x, la integral que se expresa corresponde
al área bajo la curva.

f(x)

x
a b

En la mayoría de las ocasiones los métodos numéricos para integración, al realizar un


análisis se empieza con un planteamiento desde una perspectiva gráfica.
Un planteamiento con sentido común es el de dividir el área en segmentos verticales, o
bandas, con una altura igual al valor de la función en el punto medio de cada banda.

91
Métodos numéricos

f(x)

x
a b

El área de los rectángulos se calcula y se suman para estimar el área total.


En este planteamiento, se supone que el valor de los puntos medios proporciona una
aproximación válida de la altura promedio de la función en cada banda, es posible obtener
una estimación mejor, usando más (y delgadas) bandas para aproximar la integral.
Las fórmulas de integración de NEWTON-COTES son los esquemas más comunes dentro
de la integración numérica. Se basan en la estrategia de reemplazar una función complicada
con alguna función aproximada que sea más fácil de integrar:
Las fórmulas de Newton-Cotes están conformadas por las bien conocidas reglas del
trapecio y de Simpson (regla de un tercio y de tres octavos).

Fórmulas de integración de Newton-Cotes


Estas fórmulas se basan en la idea de integrar una función polinomial en vez de f(x).
𝑏 𝑏
∫ 𝑓(𝑥)𝑑𝑥 ≈ ∫ 𝑓𝑛 (𝑥)𝑑𝑥
𝑎 𝑎

Donde fn(x)=a0 + a1x +…+ anx es un polinomio de interpolación de grado n para ciertos
n

datos de f(x) que se escogen apropiadamente.

Es importante observar que estas fórmulas se pueden aplicar inclusive a una tabla de datos,
ya que lo que se usa es un polinomio de interpolación, el cual puede ser calculado con la
tabla.

Dentro de las fórmulas de Newton-Cotes, existen las formas cerradas y abiertas. En las
formas cerradas se conocen los valores de f(a) y f(b); en caso contrario, se llaman formas
abiertas.

92
Métodos numéricos

Regla del Trapecio


Para entender el método de la regla del trapecio debe tomarse en cuenta la interpretación
geométrica de una integral como el área bajo la curva, siendo así, considérese el área de un
trapecio entre dos puntos de una curva como se muestra en la figura.

Yi-1
A1
f(a) f(b)

x
a b
(b-a)

𝒃
𝒇(𝒂) + 𝒇(𝒃)
∫ 𝒇(𝒙)𝒅𝒙 ≈ (𝒃 − 𝒂) [ ] [𝑅𝑒𝑔𝑙𝑎 𝑑𝑒𝑙 𝑇𝑟𝑎𝑝𝑒𝑐𝑖𝑜]
𝒂 𝟐

Ancho Altura promedio

A esta ecuación el nombre se debe a la interpretación geométrica que le podemos dar a la


fórmula. El polinomio de interpolación para una tabla que contiene dos datos es una línea
recta. La integral, corresponde al área bajo la línea recta en el intervalo [a, b] que es
precisamente el área del trapecio que se forma.

Ejemplo:
1 2
Aplicar la regla del Trapecio para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥
Solución:
2
Usamos la ecuación directamente con los siguientes datos: a=0, b=1 y f(x)= 𝑒 𝑥

1 2 2
𝑥2
𝑓(0) + 𝑓(1) 𝑒 (0) + 𝑒 (1) 1 + 2.7183 3.7183
∫ 𝑒 𝑑𝑥 ≈ (1 − 0) [ ]= = = = 𝟏. 𝟖𝟓𝟗𝟏𝟓
0 2 2 2 2

93
Métodos numéricos

Ejemplo:
4 𝑒𝑥
Aplicar la regla del Trapecio para aproximar la integral: ∫2 𝑑𝑥
𝑥

Solución:
Igual que en el ejemplo anterior, sustituimos los datos de manera directa en la fórmula del
𝑒𝑥
trapecio. En este caso, tenemos los datos: a=2, b=4, y 𝑓(𝑥) = 𝑥

Por lo tanto, tenemos que:


𝑒(2) 𝑒(4)
4 𝑒𝑥 𝑓(2)+𝑓(4) 2
+
4
∫2 𝑥 𝑑𝑥 ≈ (4 − 2) [ ]= = 3.6945 + 13.6495 = 𝟏𝟕. 𝟑𝟒𝟒
2 2

Algoritmo para la regla del Trapecio


function trapecio=trapecio()
clc;
f =input('Dame la función : ','s');
a =input('Dame el limite inferior: ');
b =input('Dame el limite superior: ');
ezplot(f,[a,b]);
grid,title('\bf Regla del Trapecio')
f=inline(f);
aprox=((b-a)/2)*(f(a)+f(b));
fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

Solución en Matlab
Dame la función: exp(x^2)
Dame el limite inferior : 0
Dame el limite superior: 1

Aproximación a la integral= 1.859141

94
Métodos numéricos

Ejercicios
6
1. Aplicar la regla del Trapecio para aproximar la integral: ∫1 (2 + 𝑠𝑒𝑛( 2√𝑥))𝑑𝑥
𝜋
2. Aplicar la regla del Trapecio para aproximar la integral: ∫02 𝑠𝑒𝑛2 (𝑥) 𝑑𝑥

1
3. Aplicar la regla del Trapecio para aproximar la integral: ∫0 √1 − 𝑥 3 𝑑𝑥

1 1
4. Aplicar la regla del Trapecio para aproximar la integral: ∫0 𝑑𝑥
1+𝑥 5

2 𝑒𝑥
5. Aplicar la regla del Trapecio para aproximar la integral: ∫1 𝑑𝑥
𝑥

5 2
6. Aplicar la regla del Trapecio para aproximar la integral: ∫1 𝑒 −𝑥 𝑑𝑥

2 ln 𝑥
7. Aplicar la regla del Trapecio para aproximar la integral: ∫1 𝑑𝑥
𝑥+1

𝜋
8. Aplicar la regla del Trapecio para aproximar la integral: ∫0 𝑥 𝑡𝑎𝑛𝑥 𝑑𝑥
4

3 1
9. Aplicar la regla del Trapecio para aproximar la integral: ∫2 𝑑𝑥
ln 𝑥

𝜋
10. Aplicar la regla del Trapecio para aproximar la integral: ∫02 𝑠𝑒𝑛(𝑥 2 ) 𝑑𝑥

95
Métodos numéricos

Regla del Trapecio compuesto


La regla del trapecio se puede ampliar si subdividimos el intervalo [a, b] en n subintervalos,
𝒃−𝒂
todos de la misma longitud 𝒉 = 𝑦 𝑥𝑖 = 𝑎 + 𝑖ℎ.
𝒏

Sea P= {x0, x1,…xn} la partición que se forma al hacer dicha subdivisión. Usando
propiedades de la integral tenemos que:
𝑎 𝑥1 𝑥2 𝑥𝑛
∫ 𝑓(𝑥)𝑑𝑥 = ∫ 𝑓(𝑥)𝑑𝑥 + ∫ 𝑓(𝑥)𝑑𝑥 + ⋯ + ∫ 𝑓(𝑥)𝑑𝑥
𝑏 𝑥0 𝑥1 𝑥𝑛−1

Aplicando la regla del trapecio a cada una de las integrales, obtenemos:


𝑏
𝑓(𝑥0 ) + 𝑓(𝑥1 ) 𝑓(𝑥𝑛−1 ) + 𝑓(𝑥𝑛 )
∫ 𝑓(𝑥)𝑑𝑥 ≈ (𝑥1 − 𝑥0 ) [ ] + ⋯ + (𝑥𝑛 − 𝑥𝑛−1 ) [ ]
𝑎 2 2

Sustituyendo el valor de h y usando la notación sigma, tenemos finalmente:

𝒃
𝒇(𝒙𝟎 ) + 𝟐 ∑𝒏−𝟏
𝒊=𝟏 𝒇(𝒙𝒊 ) + 𝒇(𝒙𝒏 )
∫ 𝒇(𝒙)𝒅𝒙 ≈ (𝒃 − 𝒂) [ ] [𝑇𝑟𝑎𝑝𝑒𝑐𝑖𝑜 𝑐𝑜𝑚𝑝𝑢𝑒𝑠𝑡𝑜]
𝒂 𝟐𝒏

Esta es la regla del Trapecio para n subintervalos. Obviamente, esperamos que entre más
subintervalos usemos, mejor sea la aproximación a la integral.
y

Yi-1
A1 A2 A3 A4 An

x
X0 X1

Ejemplo:
1 2
Aplicar la regla del Trapecio ompuesto para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥 en 5
subintervalos.
Solución:

96
Métodos numéricos

2
Usamos la ecuación directamente con los siguientes datos: a=0, b=1 y f(x)= 𝑒 𝑥 y en este
caso, identificamos n=5, y la partición generada es: P= {0, 0.2, 0.4, 0.6, 0.8, 1}.
Así aplicando la formula tenemos:

1
2 𝑓(0) + 2(𝑓(0.2) + 𝑓(0.4) + 𝑓(0.6) + 𝑓(08)) + 𝑓(1)
∫ 𝑒 𝑥 𝑑𝑥 ≈ (1 − 0) [ ]
0 2(5)
2 2 2 2 2 2
𝑒 (0) + 2(𝑒 (0.2) + 𝑒 (0.4) + 𝑒 (0.6) + 𝑒 (0.8) ) + 𝑒 (1)
=[ ]
10
1 + 2(1.0408 + 1.1735 + 1.4333 + 1.8965) + 2.7183 14.8065
=[ ]=
10 10
= 𝟏. 𝟒𝟖𝟎𝟔𝟓

Nota: Mencionar que el valor verdadero de esta integral es de 1.4626


Así, vemos que, con 5 intervalos, la aproximación no es tan mala, pero para un mejor
resultado se recomienda un número mayor de subintervalos.

Ejemplo:
Aplicar el método del Trapecio compuesto para aproximar:
Puntos 0 1 2 3 4 5
x 0 0.1 0.2 0.3 0.4 0.5
f(x) 1 7 4 3 5 2

Integrando por el método del Trapecio


En primer lugar, se debe verificar que el intervalo entre los puntos xi sea constante, o sea,
que el valor de “h”, sea constante.

x 0 0.1 0.2 0.3 0.4 0.5


h h h h h
(0.1) (0.1) (0.1) (0.1) (0.1)

Luego aplicando la formula, si sabemos que n=5

0.5 4

𝐼=∫ 𝑓(𝑥)𝑑𝑥 = [𝑓(𝑥) + 2 ∑ 𝑓(𝑥𝑖) + 𝑓(𝑥5)]
0 2
𝑖=0

𝑏−𝑎
ℎ= = 0.1 𝑄𝑢𝑒 𝑒𝑠 𝑒𝑙 𝑡𝑎𝑚𝑎ñ𝑜 𝑑𝑒𝑙 𝑖𝑛𝑡𝑒𝑟𝑣𝑎𝑙𝑜
𝑛

97
Métodos numéricos

Luego desarrollando la fórmula para 6 intervalos


𝐼= [𝑓(𝑥0 ) + 2[𝑓(𝑥1 ) + 𝑓(𝑥2 ) + 𝑓(𝑥3 ) + 𝑓(𝑥4 )] + 𝑓(𝑥5 )]
2

Remplazando en la fórmula:

0.1
𝐼= [1 + 2(7 + 4 + 3 + 5) + 2] = 𝟐. 𝟎𝟓
2

Algoritmo de la regla del Trapecio Compuesta


function trapeciocompuesta=trapeciocompuesta()
clc;
f =input('Dame la función : ','s');
a =input('Dame el punto a: ');
b =input('Dame el punto b: ');
n =input('Dame el número de subintervalos: ');
ezplot(f,[a,b]);
grid,title('\bf Regla del Trapecio Compuesta')
f=inline(f);

h=((b-a)/(2*n));
sumxi=0;
for i=1:n-1
x=a+h*(2*i);
sumxi=sumxi+feval(f,x);
end

aprox=((b-a)/(2*n))*(f(a)+ 2*sumxi + f(b));


fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

Solución en Matlab
Dame la función: exp(x^2)
Dame el punto a: 0
Dame el punto b: 1
Dame el número de subintervalos: 5

Aproximación a la integral= 1.480655

98
Métodos numéricos

Ejercicios
1. Aplicar la regla del Trapecio Compuesta para aproximar la integral:
2
∫0 cos (𝑥 2 )𝑑𝑥 en 3 subintervalos.

2. Aplicar la regla del Trapecio Compuesta para aproximar la integral:


𝜋
∫02 cos(𝑥 2 )𝑑𝑥 en 3 subintervalos.

3. Aplicar la regla del Trapecio Compuesta para aproximar la


13
integral: ∫0 √1 + 𝑥 3 𝑑𝑥 en 3 subintervalos.

𝜋
4. Aplicar la regla del Trapecio Compuesta para aproximar la integral:∫02 √𝑠𝑒𝑛 𝑥 𝑑𝑥
en 5 subintervalos.

5. Aplicar la regla del Trapecio Compuesta para aproximar la


13
integral:∫0 √𝑥 + 𝑥 2 𝑑𝑥 en 5 subintervalos.

1
6. Aplicar la regla del Trapecio Compuesta para aproximar la integral: ∫0 (9 −
1
𝑥 2 )3 𝑑𝑥 en 5 subintervalos.
1
7. Aplicar la regla del Trapecio Compuesta para aproximar la integral: ∫0 √tan 𝑥 𝑑𝑥
en 9 subintervalos.
1 𝑠𝑒𝑛 𝑥
8. Aplicar la regla del Trapecio Compuesta para aproximar la integral: ∫0 𝑑𝑥
𝑥
en 9 subintervalos.

9. Aplicar la regla del Trapecio Compuesta para aproximar la integral:


𝜋 1
∫0 1+𝑠𝑒𝑛2𝑥 𝑑𝑥 en 9 subintervalos.

𝜋 cos 𝑥−1
10. Aplicar l regla del Trapecio Compuesta para aproximar la integral: ∫0 𝑑𝑥
𝑥2
En 9 subintervalos.

99
Métodos numéricos

Regla de Simpson
Basado en la utilización de segmentos de parábola para aproximar los arcos de curva, en
lugar de emplear segmentos de recta; es decir utilizar curvas en lugar de una poligonal, se
obtiene una mayor precisión en el cálculo de integrales definidas.
Además de aplicar la regla trapezoidal con segmentos cada vez más finos, otra manera de
obtener una estimación más exacta de una integral es la de usar polinomios de orden
superior para conectar los puntos. Por ejemplo, si hay un punto medio extra entre f(a) y
f(b), entonces los tres puntos se pueden conectar con un polinomio de tercer orden.
A las fórmulas resultantes de calcular la integral bajo estos polinomios se les llaman Reglas
de Simpson.

Regla de Simpson 1/3


La regla de Simpson de 1/3 aproxima el área bajo la curva de f(x) mediante parábolas como
se muestra la figura, se hace pasar un polinomio de segundo orden por cada tres puntos. El
polinomio definido por los puntos xi-1, xi, y xi+1 puede obtener mediante el polinomio de
interpolación de Newton.

𝑃2 (𝑥) = 𝑎1 + 𝑎2 (𝑥 − 𝑥𝑖−1 ) + 𝑎3 (𝑥 − 𝑥𝑖−1 )(𝑥 − 𝑥𝑖−1 )

Yi-1
A1
f(a) f(b)

x
a b
xm

𝒃
𝒇(𝒂) + 𝟒𝒇(𝒙𝒎 ) + 𝒇(𝒃) 1
∫ 𝒇(𝒙)𝒅𝒙 ≈ (𝒃 − 𝒂) [ ] [𝑅𝑒𝑔𝑙𝑎 𝑑𝑒 𝑆𝑖𝑚𝑝𝑠𝑜𝑛 ]
𝒂 𝟔 3

Nota: El punto xm es el punto que está a la mitad entre a y b.

100
Métodos numéricos

Ejemplo:
1 2
Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥

Solución:
2
Aplicamos la formula directamente con los datos siguientes: a=0, b=1, xm=0.5 y f(x)= 𝑒 𝑥 .

1 2
𝑥2
𝑓(0) + 4𝑓(0.5) + 𝑓(1) 𝑒 (0) + 4𝑒 (0.5) + 𝑒 (1)
∫ 𝑒 𝑑𝑥 ≈ (1 − 0) [ ]=
0 6 6
1 + 5.1361 + 2.7183 8.8544
= = = 𝟏. 𝟒𝟕𝟓𝟕𝟑
6 6

Ejemplo:
4 𝑒𝑥
Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫2 𝑑𝑥
𝑥

Solución:
𝑒𝑥
Aplicamos la formula directamente con los datos siguientes: a=2, b=4, xm=3 y f(x)= 𝑥

𝑒2 𝑒3 𝑒4
𝑒𝑥4
𝑓(2) + 4𝑓(3) + 𝑓(4) + 4( 3) +
∫ 𝑑𝑥 ≈ (4 − 2) [ ] = [2 4
]
2 𝑥 6 3
3.6945 + 26.7807 + 13.6495 44.1247
=( )= = 𝟏𝟒. 𝟕𝟎𝟖𝟐
3 3

Algoritmo de la Regla de Simpson 1/3


function simpson13=simpson13()
clc;
f =input('Dame la función : ','s');
a =input('Dame el punto a: ');
b =input('Dame el punto b: ');
ezplot(f,[a,b]);
grid,title('\bf Regla de Simpson 1/3')
f=inline(f);

aprox=((b-a)/6)*(f(a)+4*f((a+b)/2)+f(b));
fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

101
Métodos numéricos

Ejercicios
6
1. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫1 (2 + 𝑠𝑒𝑛( 2√𝑥))𝑑𝑥
𝜋
2. Aplicar la regla Simpson 1/3 para aproximar la integral: ∫02 𝑠𝑒𝑛2 (𝑥) 𝑑𝑥

1
3. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫0 √1 − 𝑥 3 𝑑𝑥

1 1
4. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫0 𝑑𝑥
1+𝑥 5

2 𝑒𝑥
5. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫1 𝑑𝑥
𝑥

5 2
6. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫1 𝑒 −𝑥 𝑑𝑥

2 ln 𝑥
7. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫1 𝑑𝑥
𝑥+1

𝜋
8. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫0 𝑥 𝑡𝑎𝑛𝑥 𝑑𝑥
4

3 1
9. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫2 𝑑𝑥
ln 𝑥

𝜋
10. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫02 𝑠𝑒𝑛(𝑥 2 ) 𝑑𝑥

1
11. Aplicar la regla de Simpson 1/3 para aproximar la integral: ∫0 √1 + 2𝑥 𝑑𝑥

102
Métodos numéricos

Regla de Simpson 1/3 Compuesta


Al igual que con la regla del trapecio, podemos extender la regla de Simpson de 1/3, si
𝒃−𝒂
subdividimos el intervalo [a, b] en n subintervalos de la misma longitud 𝒉 = .
𝒏

Sea P={x0, x1,…,xn} la partición que se forma al hacer la subdivisión, y denotemos por xMi
ϵ [xi-1, xi] el punto medio en cada subintervalo.

Yi-1

f(a) f(b)

xm x
a b

𝒃
𝒇(𝒙𝟎 ) + 𝟒 ∑𝒏𝒊=𝟏 𝒇(𝒙𝒎𝒊 ) + 𝟐 ∑𝒏−𝟏
𝒊=𝟏 𝒇(𝒙𝒊 ) + 𝒇(𝒙𝒏 ) 1
∫ 𝒇(𝒙)𝒅𝒙 ≈ (𝒃 − 𝒂) [ ] [𝑆𝑖𝑚𝑝𝑠𝑜𝑛 𝑐𝑜𝑚𝑝𝑢𝑒𝑠𝑡𝑎]
𝒂 𝟔𝒏 3

Ejemplo:
1 2
Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥 en 5
subintervalos.

Solución:
En este caso, tenemos que n=5, y la partición que se genera es: P= {0, 0.2, 0.4, 0.6, 0.8, 1},
además los puntos medios de cada subintervalo son: PM= {0.1, 0.3, 0.5, 0.7, 0.9}
Por lo que sustituimos los datos en la ecuación para obtener:

103
Métodos numéricos

1
2
∫ 𝑒 𝑥 𝑑𝑥
0
𝑓(0) + 4[𝑓(0.1) + 𝑓(0.3) + ⋯ + 𝑓(0.9)] + 2[𝑓(0.2) + 𝑓(0.4) + ⋯ + 𝑓(0.8)] + 𝑓(1)
≈ (1 − 0) [ ]
6(5)
1 + 4[1.0101 + 1.0942 + 1.2840 + 1.6323 + 2.247] + 2[1.0408 + 1.1735 + 1.4333 + 1.8965] + 2.7183
=[ ]
30
1 + 4.0404 + 4.3768 + 5.136 + 6.5292 + 8.988 + 2.0816 + 2.347 + 2.8666 + 3.7930 + 2.7183
=[ ]
30
43.8769
= = 𝟏. 𝟒𝟔𝟐𝟔
30

Algoritmo de la regla de Simpson 1/3 compuesta


function simpson13compuesta=simpson13compuesta()
clc;
clear;
f=input('Dame la función : ','s');
a=input('Dame el punto a: ');
b=input('Dame el punto b: ');
n=input('Dame el número de subintervalos: ');

ezplot(f,[a,b]);
grid,title('\bf Regla de Simpson 1/3 compuesta')
f=inline(f);

h=(b-a)/(2*n);

sumxi=0;
for i=1:n-1
x=a+h*(2*i);
sumxi=sumxi+feval(f,x);
end

sumxmi=0;
for i=1:n
x=a+h*(2*i-1);
sumxmi=sumxmi+feval(f,x);
end

aprox=((b-a)/(6*n))*(f(a)+ 4*sumxmi + 2*sumxi + f(b));


fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

104
Métodos numéricos

Ejercicios
1. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:
2
∫0 cos (𝑥 2 )𝑑𝑥 en 3 subintervalos.

2. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


𝜋
∫02 cos(𝑥 2 )𝑑𝑥 en 3 subintervalos.

3. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la


13
integral: ∫0 √1 + 𝑥 3 𝑑𝑥 en 3 subintervalos.

4. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la


𝜋
integral:∫0 √𝑠𝑒𝑛 𝑥 𝑑𝑥 en 5 subintervalos.
2

5. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la


13
integral:∫0 √𝑥 + 𝑥 2 𝑑𝑥 en 5 subintervalos.

1
6. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral: ∫0 (9 −
1
𝑥 2 )3 𝑑𝑥 en 5 subintervalos.

7. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


1
∫0 √tan 𝑥 𝑑𝑥 en 9 subintervalos.

8. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


1 𝑠𝑒𝑛 𝑥
∫0 𝑥 𝑑𝑥 en 9 subintervalos.

9. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


𝜋 1
∫0 1+𝑠𝑒𝑛2𝑥 𝑑𝑥 en 9 subintervalos.

10. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


𝜋 cos 𝑥−1
∫0 𝑥 2 𝑑𝑥 en 9 subintervalos.

11. Aplicar la regla de Simpson 1/3 Compuesta para aproximar la integral:


2 2
∫0 𝑒 (cos 𝑥 ) 𝑑𝑥 en 9 subintervalos.

12. Aplicar la regla de Simpson 1/3 compuesta para aproximar la integral doble:
2 𝜋
∫0 ∫1 𝑑𝑥𝑑𝑦 en 9 subintervalos.

13. Aplicar la regla de Simpson 1/3 compuesta para aproximar la integral doble:
2 4
∫−2 ∫0 (𝑥 2 − 3𝑦 2 + 𝑥𝑦 3 )𝑑𝑥𝑑𝑦 en 9 subintervalos.

105
Métodos numéricos

106
Métodos numéricos

Regla de Simpson 3/8


El método de Simpson 3/8 aproxima el área bajo la curva de f(x) mediante polinomios
cúbicos. Por cada cuatro puntos se hace pasar un polinomio de tercer orden. Para puntos
xi, xi+1, xi+2, xi+3 el área bajo la curva es:

Yi-1
A1
f(a) f(b)

x
a xm xm b

𝒃
𝒇(𝒙𝟎 ) + 𝟑𝒇(𝒙𝟏 ) + 𝟑𝒇(𝒙𝟐 ) + 𝒇(𝒙𝟑 ) 3
∫ 𝒇(𝒙)𝒅𝒙 ≈ (𝒃 − 𝒂) [ ] [𝑅𝑒𝑔𝑙𝑎 𝑑𝑒 𝑆𝑖𝑚𝑝𝑠𝑜𝑛 ]
𝒂 𝟖 8

Ejemplo:
4
Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫1 𝑒 𝑥 𝐼𝑛 𝑥 𝑑𝑥

Solución:
En este caso tenemos los siguientes datos: x0=1, x1=2, x2=3, x3=4 y f(x)= 𝑒 𝑥 𝐼𝑛 𝑥, los cuales
sustituimos en la fórmula para obtener:

4
𝑓(1) + 3𝑓(2) + 3𝑓(3) + 𝑓(4)
∫ 𝑒 𝑥 𝐼𝑛 𝑥 𝑑𝑥 ≈ (4 − 1) [ ]
1 8
3
= [𝑒 1 ln(1) + 3𝑒 2 ln(2) + 3𝑒 3 ln(3) + 𝑒 4 ln(4)] = 𝟓𝟖. 𝟗𝟔𝟗𝟖
8

Nota: De x0 a x4 que son 3 segmentos, es porque se unen 4 puntos y es constante este


número de segmentos.

107
Métodos numéricos

Algoritmo de la regla de Simpson 3/8


function simpson38=simpson38()
clc;
clear;
f=input('Dame la función : ','s');
a=input('Dame el punto a: ');
b=input('Dame el punto b: ');
ezplot(f,[a,b]);
grid,title('\bf Regla de Simpson 3/8')
f=inline(f);

h=((b-a)/3);
x=a;
sum=0;
for i=2:3
x= x + h;
sum=sum + 3*f(x);
end
aprox=((b-a)/8)*(f(a)+sum +f(b));
fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

108
Métodos numéricos

Ejercicios
6
1. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫1 (2 + 𝑠𝑒𝑛( 2√𝑥))𝑑𝑥
𝜋
2. Aplicar la regla Simpson 3/8 para aproximar la integral: ∫02 𝑠𝑒𝑛2 (𝑥) 𝑑𝑥

1
3. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫0 √1 − 𝑥 3 𝑑𝑥

1 1
4. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫0 𝑑𝑥
1+𝑥 5

2 𝑒𝑥
5. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫1 𝑑𝑥
𝑥

5 2
6. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫1 𝑒 −𝑥 𝑑𝑥

2 ln 𝑥
7. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫1 𝑑𝑥
𝑥+1

𝜋
8. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫0 𝑥 𝑡𝑎𝑛𝑥 𝑑𝑥
4

3 1
9. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫2 𝑑𝑥
ln 𝑥

𝜋
10. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫02 𝑠𝑒𝑛(𝑥 2 ) 𝑑𝑥

1
11. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫0 √1 + 2𝑥 𝑑𝑥

1 2
12. Aplicar la regla de Simpson 3/8 para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥

109
Métodos numéricos

Regla de Simpson 3/8 compuesta


Al igual que en los dos casos anteriores, la regla de Simpson de 3/8, se puede extender si
𝒃−𝒂
subdividimos el intervalo [a, b] en n intervalos de la misma longitud 𝒉 = .
𝒏

Sea x0, x1, x2,…xn la partición determinada de esta forma. Cada subintervalo [xi-1, xi] lo
dividimos en 3 partes iguales, y sean yi y zi los puntos determinados así:

Yi-1

f(a) f(b)

x
a xm b

Aplicando la regla de 3/8 en cada uno de los intervalos tenemos:


𝒏 𝒏−𝟏
𝒃
𝒃−𝒂 3
∫ 𝒇(𝒙)𝒅𝒙 = [𝒇(𝒙𝟎 ) + 𝟑 (∑[𝒇(𝒚𝒊 ) + 𝒇(𝒛𝒊 )]) + 𝟐 ∑ 𝒇(𝒙𝒊 ) + 𝒇(𝒙𝒏 )] [𝑆𝑖𝑚𝑝𝑠𝑜𝑛 𝑐𝑜𝑚𝑝𝑢𝑒𝑠𝑡𝑎]
𝒂 𝟖𝒏 8
𝒊=𝟏 𝒊=𝟏

Ejemplo:
4
Aplicar la regla de Simpson 3/8 compuesta para aproximar la integral: ∫1 𝑒 𝑥 𝐼𝑛 𝑥 𝑑𝑥 en 3
subintervalos.

Solución:
Identificamos n=3 y la partición corresponde: P = {1, 2, 3, 4}

Al considerar los puntos que dividen en tres partes iguales a cada subintervalo, tenemos los
siguientes datos:

110
Métodos numéricos

1
(1 + 0.333333) = 1.333333
(1.333333 + 0.333333) = 1.666666
(1.666666 + 0.333333) = 2
(2 + 0.333333) = 2.333333
(2.333333 + 0.333333) = 2.666666
(2.666666 + 0.333333) = 3
(3 + 0.333333) = 3.333333
(3,333333 + 0.333333) = 3.666666
(3.666666 + 0.333333) = 4
𝑏−𝑎
Nota: Se utiliza la fórmula de ℎ = para sacar el valor de P y para sacar el valor de P Mi
𝑛

se utiliza 𝑝𝑚 = y aquí el 3 es constante por el número de segmentos que une los 4
3

puntos.

Sustituyendo todos los datos en la formula, obtenemos:


4
4−1
∫ 𝑒 𝑥 ln 𝑥 𝑑𝑥 ≈ [𝑓(1)
1 8(3)
+ 3[𝑓(1.333333) + 𝑓(1.666666) + 𝑓(2.333333) + 𝑓(2.666666)
3
+ 𝑓(3.333333) + 𝑓(3.666666)] + 2[𝑓(2) + 𝑓(3)] + 𝑓(4)] = [𝑒 1 ln(1)
24
+ 3[𝑒 1.333333 ln(1.333333) + ⋯ + 𝑒 3.666666 ln (3.666666)]
+ 2[𝑒 2 ln(2) + 𝑒 3 ln(3)] + 𝑒 4 ln (4)] = 𝟓𝟕. 𝟗𝟔𝟖𝟕𝟖

Algoritmo de la regla de Simpson 3/8 compuesta


function simpson38compuesta=simpson38compuesta ()
clc;
clear;
f=input('Dame la función : ','s');
a =input('Dame el punto a: ');
b =input('Dame el punto b: ');
n =input('Dame el número de subintervalos: ');

ezplot(f,[a,b]);
grid,title('\bf Regla de Simpson 3/8 compuesta')
f=inline(f);
h =(b-a)/(2*n);

f0=0;
for i=1:n-1
x=a+h*(2*i);
f0=f0+f(x);
end
f1=0;
for i=1:n

111
Métodos numéricos

x=a+h*(2*i-1);
f1=f1+f(x);
end
f0=2*f0+4*f1;
f0=f0+f(a)+f(b);
aprox=(h/3)*f0;
fprintf('\n\nAproximación a la integral=%12.6f \n',aprox);
end

112
Métodos numéricos

Ejercicios
1. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:
2
∫0 cos (𝑥 2 )𝑑𝑥 en 3 subintervalos.

2. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


𝜋
∫02 cos(𝑥 2 )𝑑𝑥 en 3 subintervalos.

3. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la


13
integral: ∫0 √1 + 𝑥 3 𝑑𝑥 en 3 subintervalos.

4. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la


𝜋
integral:∫0 √𝑠𝑒𝑛 𝑥 𝑑𝑥 en 5 subintervalos.
2

5. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la


13
integral:∫0 √𝑥 + 𝑥 2 𝑑𝑥 en 5 subintervalos.

1
6. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral: ∫0 (9 −
1
𝑥 2 )3 𝑑𝑥 en 5 subintervalos.

7. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


1
∫0 √tan 𝑥 𝑑𝑥 en 9 subintervalos

8. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


1 𝑠𝑒𝑛 𝑥
∫0 𝑥 𝑑𝑥 en 9 subintervalos.

9. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


𝜋 1
∫0 1+𝑠𝑒𝑛2𝑥 𝑑𝑥 en 9 subintervalos.

10. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


𝜋 cos 𝑥−1
∫0 𝑥 2 𝑑𝑥 en 9 subintervalos.

11. Aplicar la regla de Simpson 3/8 Compuesta para aproximar la integral:


2 2
∫0 𝑒 (cos 𝑥 ) 𝑑𝑥 en 9 subintervalos.

12. Aplicar la regla de Simpson de 3/8 Compuesta para aproximar la integral:


1 2
∫0 𝑒 𝑥 𝑑𝑥 en 9 subintervalos.

13. Aplicar la regla de Simpson 1/3 compuesta y Simpson 3/8 compuesta para
𝑥 2
aproximar la integral: ∫1 ∫1 (𝑦 + exp(𝑥 ∗ 𝑦) + 𝑥 2 )𝑑𝑥𝑑𝑦 en 10 subintervalos.

113
Métodos numéricos

14. Las áreas (A) de la sección transversal de una corriente se requieren para varias
tareas de la ingeniería de recursos hidráulicos, como el pronóstico del
escurrimiento y el diseño de presas. A menos que se disponga de dispositivos
electrónicos muy avanzados para obtener perfiles continuos del fondo del canal,
el ingeniero debe basarse en mediciones discretas de la profundidad para
calcular A.
En la figura se representa un ejemplo de sección transversal común de una
corriente. Los puntos de los datos representan ubicaciones en las que anclo un
barco y se hicieron mediciones de la profundidad. Utilice aplicaciones de la
regla del Trapecio, Simpson 1/3 y Simpson 3/8 para estimar el área de la sección
transversal representada por esos datos (h=2m).

114
Métodos numéricos

Método de Romberg
La integración de Romberg es una técnica que está especialmente diseñada para alcanzar la
eficiencia en las integrales numéricas de funciones: Está basada en la extrapolación de
Richardson, método que permite generar una estimación numérica más exacta a partir de
dos estimaciones numéricas menos exactas. En concreto, está basada en aproximaciones
obtenidas a partir de la regla trapezoidal.
𝑏
El procedimiento de Romberg para aproximar ∫𝑎 𝑓(𝑥)𝑑𝑥 consiste en lo siguiente:
Aplicando la regla del Trapecio sucesivamente para tamaños de paso hk variables, así:

𝑁𝑖𝑣𝑒𝑙 1 𝑁𝑖𝑣𝑒𝑙 2 𝑁𝑖𝑣𝑒𝑙 3 𝑁𝑖𝑣𝑒𝑙 𝑘


4 1
(𝑅𝑒𝑔𝑙𝑎 𝑑𝑒𝑙 𝑡𝑟𝑎𝑝𝑒𝑐𝑖𝑜) 𝐼(ℎ2 ) − 𝐼(ℎ1 )
3 3
𝑏−𝑎
ℎ1 = 0 ∎
2
𝑏−𝑎
ℎ2 = 1 ∎ ∎
2
𝑏−𝑎
ℎ𝑘 = 𝑘−1 ∎ ∎ ∎
2

𝟒𝒏 𝑨𝒏−𝟏 − 𝑨𝒏
𝑨𝒏 = [𝐸𝑥𝑡𝑟𝑎𝑝𝑜𝑙𝑎𝑐𝑖ó𝑛 𝑑𝑒 𝑅𝑖𝑐ℎ𝑎𝑟𝑑𝑠𝑜𝑛]
𝟒𝒏 − 𝟏

Ejemplo:
21
Aplicar el método de Romberg para aproximar la integral: ∫1 𝑥 𝑑𝑥 con h3.

A1

𝑏−𝑎 2−1
ℎ1 = = =1
20 1
2
1 1 1 1
∫ 𝑑𝑥 = ( + ) = 0.75
1 𝑥 2 1 2

A1= 0.75

A2

𝑏−𝑎 2−1
ℎ2 = = = 0.5
21 2

115
Métodos numéricos

1.5
1 0.5 1 1
∫ 𝑑𝑥 = ( + ) = 0.41667
1 𝑥 2 1 1.5
2
1 0.5 1 1
∫ 𝑑𝑥 = ( + ) = 0.29167
1.5 𝑥 2 1.5 2

A2= 0.41667 + 0.29167= 0.708333

A3

𝑏−𝑎 2−1
ℎ3 = = = 0.25
22 4
1.25
1 0.25 1 1
∫ 𝑑𝑥 = ( + ) = 0.22500
1 𝑥 2 1 1.25
1.5
1 0.25 1 1
∫ 𝑑𝑥 = ( + ) = 0.18333
1.25 𝑥 2 1.25 1.5
1.75
1 0.25 1 1
∫ 𝑑𝑥 = ( + ) = 0.15476
1.5 𝑥 2 1.5 1.75
2
1 0.25 1 1
∫ 𝑑𝑥 = ( + ) = 0.13393
1.75 𝑥 2 1.75 2

A3= 0.22500 + 0.18333 + 0.15476 + 0.13393 = 0.697024

Aplicando el Método de Romberg

4𝑛 𝐴𝑛−1 − 𝐴𝑛
𝐴𝑛 =
4𝑛 − 1
0.75
0.708333
0.697024

41 (0.708333) − 0.75 4(0.708333) − 0.75


𝐴1 = = = 0.694444
41 − 1 3

41 (0.697024) − 0.708333 4(0.697024) − 0.708333


𝐴1 = = = 0.693254
41 − 1 3

0.75
0.708333 0.694444
0.697024 0.693254

116
Métodos numéricos

42 (0.693254) − 0.694444 16(0.693254) − 0.694444


𝐴2 = = = 0.693175
42 − 1 15

H1 H2 H3
0.75
0.708333 0.694444
0.697024 0.693254 0.693175

Realizar el método de Romberg en Excel y en Matlab.

117
Métodos numéricos

Ejercicios
6
1. Aplicar el método de Romberg para aproximar la integral: ∫1 (2 + 𝑠𝑒𝑛( 2√𝑥))𝑑𝑥 en
h3.
𝜋
2. Aplicar el método de Romberg para aproximar la integral: ∫02 𝑠𝑒𝑛2 (𝑥) 𝑑𝑥 en h3.

1
3. Aplicar el método de Romberg para aproximar la integral: ∫0 √1 − 𝑥 3 𝑑𝑥 en h3.

1 1
4. Aplicar el método de Romberg para aproximar la integral: ∫0 𝑑𝑥 en h3.
1+𝑥 5

2 𝑒𝑥
5. Aplicar el método de Romberg para aproximar la integral: ∫1 𝑑𝑥 en h3.
𝑥

5 2
6. Aplicar el método de Romberg para aproximar la integral: ∫1 𝑒 −𝑥 𝑑𝑥 en h4.

2 ln 𝑥
7. Aplicar el método de Romberg para aproximar la integral: ∫1 𝑑𝑥 en h4.
𝑥+1

𝜋
8. Aplicar el método de Romberg para aproximar la integral: ∫0 𝑥 𝑡𝑎𝑛𝑥 𝑑𝑥 en h4.
4

3 1
9. Aplicar el método de Romberg para aproximar la integral: ∫2 𝑑𝑥 en h4.
ln 𝑥

𝜋
10. Aplicar el método de Romberg para aproximar la integral: ∫02 𝑠𝑒𝑛(𝑥 2 ) 𝑑𝑥 en h5.

1
11. Aplicar el método de Romberg para aproximar la integral: ∫0 √1 + 2𝑥 𝑑𝑥 en h5.

1 2
12. Aplicar el método de Romberg para aproximar la integral: ∫0 𝑒 𝑥 𝑑𝑥 en h5.

118
Métodos numéricos

SOLUCIÓN NUMÉRICA DE ECUACIONES DIFERENCIALES ORDINARIAS


Las leyes que gobiernan los fenómenos de la naturaleza se expresan habitualmente en
forma de ecuaciones diferenciales. Las ecuaciones del movimiento de los cuerpos (Segunda
Ley de Newton) es una ecuación diferencial de segundo orden, como lo es la ecuación que
describe los sistemas oscilantes, la propagación de las ondas, la transmisión del calor, la
difusión, el movimiento de partículas subatómicas, entre otras.

Pocas ecuaciones diferenciales tienen una solución analítica sencilla, la mayor parte de las
veces es necesario realizar aproximaciones, estudiar el comportamiento del sistema bajo
ciertas condiciones. Así, en un sistema tan simple como un péndulo, la amplitud de la
oscilación ha de ser pequeña y el rozamiento ha de ser despreciable, para obtener una
solución sencilla que describa aproximadamente su movimiento periódico.

Existen una gran diversidad de métodos numéricos para la resolución de un Problema de


Valor Inicial (PVI), con distintas características. Estos métodos se agrupan en dos
familias: los métodos de un paso y los métodos lineales multipaso.

Métodos de un paso: Se caracterizan porque el valor aproximado 𝑦𝑖 de la solución en el


punto xi se obtiene a partir del valor 𝑦𝑖−1 obtenido anteriormente.

Métodos lineales multipaso: Utilizan para el cálculo del valor aproximado yi no solo el
valor de yi-1 obtenido anteriormente, sino también los valores yi-1,….yi-j obtenidos en etapas
previas.

Una ecuación diferencial es una ecuación en donde aparecen funciones, sus derivadas,
una o más variables independientes y una o más variables dependientes, sin embargo
“ecuaciones en derivadas” sería más descriptivo. Estas se dividen en dos grupos:
1. Ecuaciones Diferenciales Ordinarias (EDO). En donde aparece sólo una variable
independiente (que se denota con x).
2. Ecuaciones Diferenciales Parciales (EDP). En las que aparece más de una variable
independiente.

119
Métodos numéricos

Método de Euler-Cauchy
El método de Euler, llamado así en honor de Leonhard Euler, es un procedimiento de
integración numérica para resolver ecuaciones diferenciales ordinarias a partir de un valor
inicial dado.

El método de Euler es el más simple de los métodos numéricos para resolver un problema
del siguiente tipo, también llamados problemas de valor inicial:
𝑦′ = 𝑓(𝑥, 𝑦)
Problema de valor inicial ={ 𝑦(𝑥0 ) = 𝑦0
𝑦(𝑥𝑖 ) =?
El método consiste en multiplicar los intervalos que van de x0 a xi en n subintervalos de
𝑥𝑖 −𝑥0
ancho h, o sea: ℎ =
𝑛

De manera que se obtiene un conjunto discreto de 𝑛 + 1 puntos: 𝑥0 , 𝑥1 , 𝑥2 , … . . 𝑥𝑛 del


intervalo de interés [𝑥0 , 𝑥𝑓 ], para cualquiera de estos puntos se cumple que:
𝑥𝑖 = 𝑥0 + 𝑖ℎ, 0 ≤ 𝑖 ≤ 𝑛
La condición inicial 𝑦(𝑥0 ) = 𝑦0 , representa el punto 𝑃0 = (𝑥0 , 𝑦0 ) por donde pasa la curva
solución de la ecuación del planteamiento inicial, la cual se denotará 𝑓(𝑥) = 𝑦.

Ya teniendo el punto 𝑃0 se puede evaluar la primera derivada de 𝑓(𝑥) en ese punto, por lo
𝑑𝑦
tanto: 𝑓 ′ (𝑥) = 𝑑𝑥 | 𝑃0 = 𝑓(𝑥0 , 𝑦0 )

y Curva de solución

(x1, y(x1))
error

(x1, y1)
Recta tangente = y’0

(x0, y0)
x0 x1 = x 0 + h x

h
Figura 8. Esquema del método de Euler

120
Métodos numéricos

Curva de solución
y (x3, y3)
(x2, y2)
(x4, y4)
(x1, y1) error

(x0, y0)

x0 x1 x2 x3 x4 x

Se resuelve para y1:


𝑦1 = 𝑦0 + (𝑥1 − 𝑥0 )𝑓(𝑥0 , 𝑦0 ) = 𝑦0 + ℎ𝑓(𝑥0 , 𝑦0 )

Es evidente que la ordenada 𝑦1 calculada de esta manera no es igual a 𝑓(𝑥1 ), pues existe un
pequeño error. Sin embargo, el valor 𝑦1 sirve para que se aproxime 𝑓’(𝑥) en el punto 𝑃 =
(𝑥1 , 𝑦1 ) y repetir el procedimiento anterior a fin de generar la sucesión de aproximaciones
siguiente:
𝑦1 = 𝑦0 + ℎ𝑓(𝑥0 , 𝑦0 )
𝑦2 = 𝑦1 + ℎ𝑓(𝑥1 , 𝑦1 )
.
.
.
𝒚𝒊+𝟏 = 𝒚𝒊 + 𝒉𝒇(𝒙𝒊 , 𝒚𝒊 )

𝒙𝒊+𝟏 ≈ 𝒙𝒊 + 𝒉
Ejemplo:
Aplicar el método de Euler para estimar un valor aproximado de la solución en 𝑥 = 1.5 de
𝑦’ = 𝑥 + 𝑦 2 , 𝑦(1) = 0 y el tamaño de paso 0.1.

En este caso se tiene 𝑦(1) = 0, por lo que el valor de 𝑥0 = 1 y 𝑦0 = 0; calcular n,


considerando el tamaño de paso de 0.1.
Despejando n de:
𝑥𝑖 − 𝑥0
ℎ=
𝑛
quedaría como:

121
Métodos numéricos

𝑥𝑖 − 𝑥0
𝑛=( )

1.5 − 1
𝑛=( )=5
0.1

𝑥0 = 𝟏
𝑦0 = 𝟎
𝑓(𝑥0 , 𝑦0 ) = (𝑥 + 𝑦 2 ) = (1 + 02 ) = 1

𝑥1 = 𝑥0 + ℎ
𝑥1 = 1 + 0.1 = 1.1
𝑦1 = 𝑦0 + ℎ𝑓(𝑥0 , 𝑦0 )
𝑦1 = 0 + 0.1 ∗ (1) = 0.1
𝑓(𝑥1 , 𝑦1 ) = (𝑥 + 𝑦 2 ) = (1.1 + 0.12 ) = 1.11

𝑥2 = 𝑥1 + ℎ
𝑥2 = 1.1 + 0.1 = 1.2
𝑦2 = 𝑦1 + ℎ𝑓(𝑥1 , 𝑦1 )
𝑦2 = 0.1 + 0.1 ∗ (1.11) = 0.211
𝑓(𝑥2 , 𝑦2 ) = (𝑥 + 𝑦 2 ) = (1.2 + 0.2112 ) = 1.244521

Así sucesivamente hasta que n sea igual a 5

n xi yi f(xi,yi)
0 1 0 1
1 1.1 0.1 1.11
2 1.2 0.211 1.244521
3 1.3 0.3354521 1.41252811
4 1.4 0.47670491 1.62724757
5 1.5 0.63942967 1.9088703

yi
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0.9 1 1.1 1.2 1.3 1.4 1.5 1.6

122
Métodos numéricos

Ejercicios

1. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =


𝑠𝑖𝑛𝑥 − 𝑙𝑛𝑦, 𝑦(0.13) = 0.32, 𝑦(0.14) =?, calcular el tamaño de paso h tomando en
cuenta que el valor de divisiones es de 4.

2. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =


−2𝑥 3 + 12𝑥 2 − 20𝑥 + 8.5, con 𝑥0 = 0, 𝑥 = 4, 𝑦0 = 1, tamaño de paso h de 0.5.

3. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =


1
𝑥 + 5𝑦 con 𝑥0 = 0, 𝑥 = 5, 𝑦0 = −3, calcular n considerando el tamaño de paso de 1.

4. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =


𝑦
√2𝑥+1 con 𝑥0 = 0, 𝑥 = 5, 𝑦0 = 4, calcular n considerando el tamaño de paso de

0.5.
5. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =
2𝑦 − 2𝑥 − 1 con h=1, y (0) =3; y (1) =?
6. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦 ′ =
3𝑒 −𝑥 − 0.5𝑦 con h=1.5, y (0) =6; y (3) =?
7. Aplicar el método de Euler para estimar un valor aproximado de la solución de 𝑦´ =
𝑥𝑦 + √𝑥 con y (0) =1, y (0.5) =?

123
Métodos numéricos

Método de Runge-Kutta (Cuarto orden)


El desarrollo de estas fórmulas se inició con el trabajo de Carl Runge en 1895 y lo
continuo M. Kutta en 1901. El método Runge-Kutta de orden 4 es la forma de los métodos
de Runge-Kutta de uso más común y así mismo más exactos para obtener soluciones
aproximadas de ecuaciones diferenciales.

El método de Runge-Kutta de cuarto orden consiste en determinar constantes apropiadas de


modo que una fórmula como:
𝑦𝑖+1 = 𝑦𝑖 + 𝑎𝑘1 + 𝑏𝑘2 + 𝑐𝑘3 + 𝑑𝑘4
Coincide con un desarrollo de Taylor hasta el término h4, es decir, hasta el quinto término.

Las ecuaciones del método de Runge-Kutta de cuarto orden son las siguientes:


𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 2𝑘2 + 2𝑘3 + 𝑘4 )
6

Donde:

𝑘1 = 𝑓(𝑥𝑖 + 𝑦𝑖 )

ℎ ℎ
𝑘2 = 𝑓 (𝑥𝑖 + , 𝑦𝑖 + 𝑘1 )
2 2
ℎ ℎ
𝑘3 = 𝑓 (𝑥𝑖 + , 𝑦𝑖 + 𝑘2 )
2 2

𝑘4 = 𝑓(𝑥𝑖 + ℎ, 𝑦𝑖 + ℎ𝑘3 )

Ejemplo:

Resolver por el método de Runge-Kutta de cuarto orden el problema de valor inicial

𝑦´ = 𝑥 2 − 3𝑦; 𝑦(0) = 1 en el intervalo 0 ≤ 𝑥 ≤ 0.4 𝑐𝑜𝑛 ℎ = 0.1.

𝑥0 = 𝟎
𝑦0 = 𝟏
𝑓(𝑥, 𝑦) = 𝑥 2 − 3𝑦

𝑘1 = 𝑓(𝑥0 , 𝑦0 ) = 𝑥 2 − 3𝑦 = 02 − 3 ∗ 1 = −3

ℎ ℎ 0.1 2 0.1
𝑘2 = 𝑓 (𝑥0 + , 𝑦0 + 𝑘1 ) = (0 + ) − 3 ∗ (1 + ( ) ∗ −3) = −2.5475
2 2 2 2

124
Métodos numéricos

ℎ ℎ 0.1 2 0.1
𝑘3 = 𝑓 (𝑥0 + , 𝑦0 + 𝑘2 ) = (0 + ) − 3 ∗ (1 + ( ) ∗ −2.5475) = −2.615375
2 2 2 2

𝑘4 = 𝑓(𝑥0 + ℎ, 𝑦0 + ℎ𝑘3 ) = (0 + 0.1)2 − 3 ∗ (1 + (0.1 ∗ −2.615375)) = −2.2053875


𝑦1 = 𝑦0 + (𝑘1 + 2𝑘2 + 2𝑘3 + 𝑘4 )
6
0.1
=1+ (−3 + 2 ∗ (−2.5475) + 2 ∗ (−2.615375) − 2.2053875)
6
= 0.741147708

Así hasta el 4 termino:

n x k1 k2 k3 k4 y
0 0 1
1 0.1 -3 -2.5475 -2.615375 -2.2053875 0.74114771
2 0.2 -2.21344313 -1.86892666 -1.92060413 -1.60726189 0.5511516
3 0.3 -1.6134548 -1.34893658 -1.38861431 -1.1468705 0.41389448
4 0.4 -1.15168344 -0.94643093 -0.9772188 -0.7885178 0.31743614

Y
1.2
1
0.8
0.6
0.4
0.2
0
0 0.1 0.2 0.3 0.4 0.5

Ejemplo:

Resolver por el método de Runge-Kutta de cuarto orden el problema de valor inicial con
una aproximación del valor de 𝑦(1.5), tomando un paso ℎ = 0.1.
𝑦 ′ = 2𝑥𝑦
{
𝑦(1) = 1

Despejando n de:

125
Métodos numéricos

𝑥𝑖 − 𝑥0
ℎ=
𝑛
quedaría como:
𝑥𝑖 − 𝑥0
𝑛=( )

1.5 − 1
𝑛=( )=5
0.1
𝑥0 = 𝟏
𝑦0 = 𝟏
𝑥𝑖 = 𝟏. 𝟓
𝑓(𝑥, 𝑦) = 2𝑥𝑦

𝑘1 = 𝑓(𝑥0 , 𝑦0 ) = 2𝑥𝑦 = 2 ∗ 1 ∗ 1 = 2

ℎ ℎ 0.1 0.1
𝑘2 = 𝑓 (𝑥0 + , 𝑦0 + 𝑘1 ) = 2 ∗ (1 + ) ∗ (1 + ( ) ∗ 2) = 2.31
2 2 2 2
ℎ ℎ 0.1 0.1
𝑘3 = 𝑓 (𝑥0 + 2 , 𝑦0 + 2 𝑘2 ) = 2 ∗ (1 + 2
)∗ (1 + ( 2 ) ∗ 2.31) = 2.34255

𝑘4 = 𝑓(𝑥0 + ℎ, 𝑦0 + ℎ𝑘3 ) = 2 ∗ (1 + 0.1) ∗ (1 + (0.1 ∗ 2.34255)) = 2.715361

ℎ 0.1
𝑦1 = 𝑦0 + (𝑘1 + 2𝑘2 + 2𝑘3 + 𝑘4 ) = 1 + (2 + 2 ∗ 2.31 + 2 ∗ 2.34255 + 2.715361)
6 6
= 1.23367435

Así hasta el 5 termino:

n x k1 k2 k3 k4 Y
0 1 1
1 1.1 2 2.31 2.34255 2.715361 1.23367435
2 1.2 2.71408357 3.14957062 3.19965163 3.72873483 1.5526954
3 1.3 3.72646896 4.34754711 4.42518188 5.18755532 1.99368677
4 1.4 5.1835856 6.08273833 6.20412395 7.31947766 2.61163323
5 1.5 7.31257305 8.63405947 8.825675 10.4826022 3.49021064

126
Métodos numéricos

4
3.5
3
2.5
2
1.5
1
0.5
0
0.9 1 1.1 1.2 1.3 1.4 1.5 1.6

Ejercicios
1. Resolver por el método de Runge-Kutta de cuarto orden el problema de valor inicial
𝑦
𝑦′ = 2 −
{ 𝑥
𝑦(1) = 0

en el intervalo 1 ≤ 𝑥 ≤ 2.5 𝑐𝑜𝑛 ℎ = 0.1

2. Resolver por el método de Runge-Kutta de cuarto orden el problema de valor inicial


con una aproximación del valor de 𝑦(1.5) para el siguiente problema de valor
inicial, tomando un paso ℎ = 0.1
𝑦 ′ = (𝑦 + 1)(𝑥 + 1) cos (𝑥 2 + 2𝑥)
{
𝑦(0) = 4

3. Resolver por el método de Runge-Kutta de cuarto orden el problema de valor inicial


𝑦´ = 𝑦𝑥 2 − 1.2𝑦; 𝑦(0) = 1

en el intervalo 0 ≤ 𝑥 ≤ 2 𝑐𝑜𝑛 ℎ = 0.5.

4. La ecuación diferencial que rige la vibración de un amortiguador en un auto de


𝑑2 𝑥 𝑑𝑥
acuerdo con la Ley de Hooke está dada por la expresión: 𝑚 𝑑2 𝑡 + 𝑐 𝑑𝑡 + 𝑘𝑥 = 0

Donde m es la masa del auto (en kg); c es un coeficiente de amortiguamiento


viscoso (Ns/m) y k es la constante de amortiguamiento del resorte en (N/m).
Transforme la expresión en una EDO, para valores de t= [0, 0.4] con tamaño de

127
Métodos numéricos

paso h=0.1 donde m=1.2x106, c=1x107 e k=1.5x109 en el caso que x=0.5 y dx/dt=0
en t=0. La ecuación analítica es la siguiente:
𝑥(𝑡) = 𝑒 −4.16667𝑡 (0.0593391 sin(35.109𝑡) + 0.5 cos(35.109𝑡))

Determinar la solución numérica x(t) utilizando el método de Runge Kutta 4 orden

128
Métodos numéricos

REFERENCIAS
1. Burden R.L. y Faires J.D. (2017), Análisis numérico, México: International Thomson
Editores
2. Chapra, S. y Canale, R. (2015), Métodos numéricos para ingenieros, México, McGraw
Hill.
3. Gerald, C. y Wheatley, P. (2004), Applied Numerical Analysis, USA: Pearson.
4. Gutiérrez, J.A., Olmos, M.A. y Casillas J.M. (2010), Análisis numérico, México:
McGraw Hill/Interamericana Editores
5. Kolman, B. y Hill, D.R. (2006), Algebra lineal, México: Pearson Educación
6. Kiusalaas, J. (2010), Numerical Methods in Engineering with MATLAB, USA:
Cambridge University Press
7. Nieves, A. y Domínguez, F.C. (2014), Métodos numéricos aplicados a la ingeniería,
México: Grupo Editorial Patria.

129

También podría gustarte