Trabajo de Metodos

INSTITUTO POLITECNICO NACIONAL
Escuela Superior de Ingeniera Mecnica y

Elctrica
Curso Propedutico para la Maestra en Ciencias en
Ingeniera Elctrica
[PROGRAMACIN Y MTODOS NUMRICOS]

PROFESOR: M.C. JESUS REYES GARCIA
ALUMNO: CARLOS ALFREDO HARRISON ROJAS
TRABAJO:
RECOPILACION DE LOS PROGRAMAS VISTOS EN CLASES
FECHA 14/DICIEMBRE/2015
Contenido
INTRODUCCIN.......................................................................................................... 2
1.
RAICES DE ECUACIONES....................................................................................... 3
1.1.
1.1.1.
EL MTODO DE BISECCIN.........................................................................3
1.1.1.1.
PROGRAMA DE TABULACION Y DE BISECCION DE UN POLINOMIO........4
1.1.1.2.
CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN....................7
1.1.2.
METODO DE REGLA FALSA........................................................................11
1.1.2.1.
PROGRAMA DE TABULACIN Y DEL MTODO DE REGLA FALSA..........12
1.1.2.2.
CAPTURA DE PANTALLA DEL PROGRAMA EN EJECUCIN....................15
1.2.
1.2.1.
MTODOS ABIERTOS.................................................................................... 17
MTODO DE NEWTON RAPHSON............................................................17
1.2.1.1.
PROGRAMA Y DIAGRAMA DE FLUJO.....................................................18
1.2.1.2.
CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN..................21
1.2.2.
MTODO DE LA SECANTE.........................................................................22
1.2.2.1.
1.2.2.2.
1.3.
1.3.1.
2.
Mtodos cerrados.......................................................................................... 3
RACES DE POLINOMIOS............................................................................... 27
Mtodo de Lin Bairstow............................................................................ 27
1.3.1.1.
1.3.1.2.
SISTEMAS DE ECUACIONES ALGEBRAICAS LINEALES.........................................33

2.1.
Eliminacin Gaussiana Simple y Sustitucin hacia atrs..............................33
INTRODUCCIN
Los mtodos numricos son tcnicas mediante las cuales es posible formular problemas de tal
forma que puedan resolverse usando operaciones aritmticas. (Steven Chapra, 2014)
Aunque hay muchos tipos de mtodos numricos, todos comparten una caracterstica comn:
Invariablemente los mtodos numricos llevan a cabo un buen nmero de clculos aritmticos.
No es raro que con el desarrollo de computadoras digitales eficientes y rpidas, el papel de los
mtodos numricos en la solucin de problemas de ingeniera haya aumentado
considerablemente en los ltimos aos.
Los mtodos numricos son herramientas extremadamente poderosas para la solucin de
problemas. Son capaces de manejar sistemas de ecuaciones grandes, no linealidades y
geometras complicadas que son comunes en la prctica de la ingeniera y que, a menudo, son
imposibles de resolver analticamente. Por lo tanto, amplan la habilidad de quien los estudia
para resolver problemas.Hay muchos problemas que no pueden plantearse al emplear
programas hechos. Si se est versado en los mtodos numricos y se es un adepto de la
programacin de computadoras, entonces se tiene la capacidad de disear programas propios
para resolver los problemas, sin tener que comprar un software costoso. En el presente
documento se usar la programacin en lenguaje C para construir nuestros propios programas
para la resolucin de problemas usando los mtodos numricos que se describen a
continuacin:
1. Races de ecuaciones. Estos problemas estn relacionados con el valor de una variable
o de un parmetro que satisface una ecuacin ya sea trascendental o un polinomio
2. Sistemas de ecuaciones algebraicas lineales. En esencia, estos problemas son
similares a los de races de ecuaciones en el sentido de que estn relacionados con
valores que satisfacen ecuaciones. Sin embargo, a diferencia de satisfacer una sola
ecuacin, se busca un conjunto de valores que satisfaga simultneamente a un
conjunto de ecuaciones algebraicas. Las ecuaciones lineales simultneas surgen en el
contexto de una variedad de problemas y en todas las disciplinas de la ingeniera.
3. Ajuste de curvas. Con frecuencia se presentar la oportunidad de ajustar curvas a un
conjunto de datos representados por puntos. Las tcnicas que se han desarrollado para
este fin pueden dividirse en dos categoras generales: regresin e interpolacin. La
regresin se emplea cuando hay un grado significativo de error asociado a los datos;
frecuentemente los resultados experimentales son de esta clase. Para estas
situaciones, la estrategia es encontrar una curva que represente la tendencia general de
los datos sin necesidad de tocar los puntos individuales. En contraste, la interpolacin
se maneja cuando el objetivo es determinar valores intermedios entre datos que estn
relativamente libres de error.
4. Ecuaciones diferenciales ordinarias. Las ecuaciones diferenciales ordinarias tienen un
enorme significado en la prctica de la ingeniera. Esto se debe a que muchas leyes
fsicas estn expresadas en trminos de la razn de cambio de una cantidad ms que
en trminos de su magnitud. Entre los ejemplos se observan desde los modelos de
prediccin demogrfica (razn de cambio de la poblacin) hasta la aceleracin de un
cuerpo en descenso (razn de cambio de la velocidad).
1. RAICES DE ECUACIONES
1.1.
Mtodos cerrados.
Los mtodos cerrados consisten en dar 2 puntos tales que stos encierren el valor de la raz,
dichos puntos es cuando al sustituirlos en la funcin, sta cambia de signo. En ste caso se
aplica un anlisis de la curva mediante la implementacin de un programa tabula.cpp en c++, el
cual evala la funcin para valores obtenidos a travs de un intervalo dado y un nmero de
puntos en los cuales se divide dicho intervalo. El programa tabula.cpp muestra al final los
puntos para los cuales en la funcin hay un cambio de signo, estos puntos se usan en los
mtodos cerrados de biseccin y regla falsa.
1.1.1.
EL MTODO DE BISECCIN.
Se sabe que un polinomio de grado n tiene n races las cuales pueden ser:
Reales y distintas.
Reales e iguales.
Complejas conjugadas.
De acuerdo a ello un polinomio de grado impar tendr por lo menos una raz real, para dichas
races ocuparemos los mtodos numricos de Biseccin, Regla Falsa, Newton Raphson.
En general, si f(x) es real y contina en el intervalo que va desde x1 hasta x2 y f(x1) y f(x2)
tienen signos opuestos, es decir,
f (x1) < 0 < f(x2)
Entonces hay al menos una raz real entre x1 y x2. Los mtodos de bsqueda incremental
aprovechan esta caracterstica localizando un intervalo en el que la funcin cambie de signo.
Entonces, la localizacin del cambio de signo (y, en consecuencia, de la raz) se logra con ms
exactitud al dividir el intervalo en varios subintervalos. Se investiga cada uno de estos
subintervalos para encontrar el cambio de signo. El proceso se repite y la aproximacin a la
raz mejora cada vez ms en la medida que los subintervalos se dividen en intervalos cada vez
ms pequeos. Por lo que sabemos que existe, al menos, una raz real. A partir de este punto
se va reduciendo el intervalo sucesivamente hasta hacerlo tan pequeo como exija la precisin
que hayamos decidido emplear. En seguida se explicara a fondo y paso a paso el mtodo:
Figura 1. Mtodo de biseccin

Criterio de convergencia.
Se debe desarrollar un criterio objetivo para decidir cundo debe terminar el mtodo. Una
sugerencia inicial sera finalizar el clculo cuando el error verdadero se encuentre por debajo
de algn nivel prefijado. Puede decidirse que el mtodo termina cuando se alcance un error
ms bajo, por ejemplo, al 0.1%. Se requiere estimar el error de forma tal que no se necesite el
conocimiento previo de la raz, se puede calcular el error relativo porcentual ea de la siguiente
manera:
Donde xr nuevo es la raz en la iteracin actual y xr anterior es el valor de la raz en la iteracin

anterior. Se utiliza el valor absoluto, ya que por lo general importa slo la magnitud de ea sin
considerar su signo. Cuando ea es menor que un valor previamente fijado es, termina el
clculo.
1.1.1.1. PROGRAMA DE TABULACION Y DE BISECCION DE UN POLINOMIO

A continuacin se presenta el cdigo del programa el cual al principio pide los
intervalos en los cuales se encuentra la raz tabulando los resultados y facilitando la
aproximacin a la raz. El ejemplo del programa se har con la siguiente funcin.
funcion=4.0 x 330 .0 x2 +70.0 x50.0
PROGRAMA
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
using namespace std;
float f(float);
void escvec2(int , float[50],char [5], float [50], char [5]);
void Tabula(float, float, int, float [50], float [50]);
int bis(float ,float , float , int, float*);
int main(int argc, char *argv[])
{
int nit,con,i, np, nr;
float raiz,x1,x2,Eps,res,xa, xb,x[50],y[50];
cout<<"Metodo de la Biseccion \n"<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de la Biseccion \n"<< endl;
cout<<"
Tabular el Polinomio \n"<< endl;
cout<<"xa: Intervalo inicial xb: Intervalo final np: Numero de partes \n"
<<endl;
cout<<"xa= ";cin>>xa;
cout<<"xb= ";cin>>xb;
cout<<"np= ";cin>>np;
Tabula(xa, xb,np,x,y);
escvec2(np+1,x,"x",y,"y");
system("PAUSE");
system("cls");
cout<<"
Metodo de la Biseccion \n"<< endl;
cout<<"Encontrar las raices \n"<< endl;
cout<<"Numero de raices del Polinomio= ";cin>>nr;
system("PAUSE");
system("cls");
for(i=1;i<=nr;i++)
{
cout<<"Raiz"<<"["<<i<<"]"<<endl;
cout<<"x1= ";cin>>x1;
cout<<"Eps= ";cin>>Eps;
cout<<"nit= ";cin>>nit;
if (bis (x1,x2,Eps,nit,&raiz)==1)
{
cout<<"\n \n Para el intervalo x1: "<<x1<<"
Para el intervalo x2:
"<<x2<<endl;
cout<<"\n El proceso converge: "<<raiz<<"\n \n \n"<<endl;
}
else
{
"<<x2<<endl;
cout<<"\n \n El proceso no converge"<< endl;
system("PAUSE");
system("cls");
}
system("PAUSE");
system("cls");
}
}
int bis(float x1,float x2, float Eps, int nit, float *res)
{
int con,i;
float x3,x4,y1,y2,y3, Ea;
con=0;
printf ("\n nit x1
x2
x3 |Xi+1 - Xi| Ea
y1
y2
y3 ");
for(i=1;i<=nit;i++)
{
x3=(x1+x2)/2.0;
if(f(x3)==0.0)
{
*res=x3;
con=1;
break;
}
y1=f(x1);
y2=f(x2);
y3=f(x3);
if(i==1)
printf("\t \n %2d %8.4f %8.4f %8.4f \t \t \t %8.4f %8.4f
%8.4f",i,x1,x2,x3,y1,y2,y3);
else
{
Ea=fabs((x3-x4)/x3)*100;
printf("\t \n \n %2d %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f",
i,x1,x2,x3,fabs(x3-x4),Ea,y1,y2,y3);
if(fabs(x3-x4)<=Eps)
{
*res=x3;
con=1;
break;
}
}
x4=x3;
if(f(x3)*f(x1)>0)
x1=x3;
else
x2=x3;
}
if(con==1)
return (con);
system("PAUSE");
}
void Tabula(float xa, float xb, int np, float x[50], float y[50])
{
int i;
float x1, dx;
x1 = xa;
dx=(xb-xa)/np;
for (i=1; i<=np ;i++)
{
x[i]=x1;
y[i]= f(x[i]);
x1=x1+dx;
}
}
void escvec2(int n, float x[50],char nomx[5], float y[50], char nomy[5])
{
int i;
for(i=1; i<=n; i++)
{
cout << nomx << "[" <<i<<"]="<<x[i]<<" "<<nomy<< "["
<<i<<"]="<<y[i]<< endl;
}
}
float f(float x)
{
return (4.0*pow(x,3)-30.0*pow(x,2)+70*x-50);
}
1.1.1.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN
Figura 2. El programa solicita intervalos y en cuantas partes hace la tabulacin.
Figura 3. Ventana de resultados de la tabulacin del polinomio.
Figura 4. Ya con la funcin tabulada podemos encontrar un intervalo menor para

encontrar la raz ahora el programa te pide los intervalos el nmero de
iteraciones y tambin el psilon.
Figura 5. Resultados del programa el mtodo converge en 16 iteraciones y nos

muestra la raz. Que para el caso de esta funcin es 1.3819.
1.1.2.
METODO DE REGLA FALSA
La falsa posicin es una alternativa basada en una visualizacin grfica. Un inconveniente del
mtodo de biseccin es que al dividir el intervalo de xl a xu en mitades iguales, no se toman en
consideracin las magnitudes de f(xl) y f(xu). Por ejemplo, si f(xl) est mucho ms cercana a
cero que f(xu), es lgico que la raz se encuentre ms cerca de xl que de xu. Un mtodo
alternativo que aprovecha esta visualizacin grfica consiste en unir f(xl) y f(xu) con una lnea
recta. La interseccin de esta lnea con el eje de las x representa una mejor aproximacin de la
raz. El hecho de que se reemplace la curva por una lnea recta da una falsa posicin de la
raz; de aqu el nombre de mtodo de la falsa posicin, o en latn, regula falsi. Tambin se le
conoce como mtodo de interpolacin lineal.Usando tringulos semejantes, la interseccin de
la lnea recta con el eje de las x se estima mediante.
al despejar Xr
Figura 6. Formula de la regla falsa
sta es la frmula de la falsa posicin. El valor de xr calculado con la ecuacin, reemplazar,
despus, a cualquiera de los dos valores iniciales, xl o xu, y da un valor de la funcin con el
mismo signo de f(xr). De esta manera, los valores xl y xu siempre encierran la verdadera raz.
El proceso se repite hasta que la aproximacin a la raz sea adecuada. El algoritmo es idntico
al de la biseccin, excepto en que la ecuacin se usa en el paso 2. Adems, se usa el mismo
criterio de terminacin para concluir los clculos.
Figura 7. Representacin grfica del mtodo de la falsa posicin. Con los tringulos
semejantes sombreados se obtiene la frmula para el mtodo.
1.1.2.1. PROGRAMA DE TABULACIN Y DEL MTODO DE REGLA FALSA

Con el fin de hacer prctico y corroborar el resultado arrojado por el mtodo de
biseccin se utilizara los mismos intervalos y la misma funcin.
funcion=4.0 x 330 .0 x2 +70.0 x50.0

PROGRAMA
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
float f1(float);
int regfal(float ,float , float , int, float*);
{
cout<<"Metodo de la Regla Falsa \n"<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de la Regla Falsa \n"<< endl;
cout<<"
<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de la Regla Falsa \n"<< endl;
system("PAUSE");
system("cls");
for(i=1;i<=nr;i++)
{
if (regfal (x1,x2,Eps,nit,&raiz)==1)
{
"<<x2<<endl;
}
else
{
"<<x2<<endl;
system("PAUSE");
system("cls");
}
system("PAUSE");
system("cls");
}
}
int regfal(float xa, float xb, float Eps, int nit, float *res)
{
int i, con;
float x1, x2, x3, x4, Ea, y1, y2, y3;
con=0;
x1=xa;
x2=xb;
x4=0;
printf(" It
x1
x2
x3 |xi+1 - xi| Ea
y1
y2
y3");
for(i=1;i<=nit;i++)
{
y1=f1(x1);
y2=f1(x2);
x3= x2 - ( ((f1(x2))*(x1-x2)) / ((f1(x1))-(f1(x2))) );
y3=f1(x3);
if(y3==0)
{
*res=x3;
con=1;
break;
}
if(i==1)
{
printf("\n %2d %8.4f %8.4f %8.4f \t \t \t %8.4f %8.4f %8.4f", i ,
x1, x2, x3, y1, y2, y3);
}
else
{
Ea=fabs((x3-x4)/x3)*100;
printf("\n %2d %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f", i ,
x1, x2, x3, fabs(x2-x3), Ea,y1, y2, y3);
{
*res=x3;
con=1;
break;
}
if((y1)*(y3)>0.0)
{
x1=x3;
}
else
{
x2=x3;
}
x4=x3;
}
}
return(con);
}
{
int i;
float x1, dx;
x1 = xa;
dx=(xb-xa)/np;
for (i=1; i<=np ;i++)
{
x[i]=x1;
y[i]= f1(x[i]);
x1=x1+dx;
}
}
{
int i;
for(i=1; i<=n; i++)
{
cout << nomx << "[" <<i<<"]="<<x[i]<<"
"<<nomy<< "["
<<i<<"]="<<y[i]<< endl;
}
}
float f1(float x)
{
}
1.1.2.2. CAPTURA DE PANTALLA DEL PROGRAMA EN EJECUCIN
Figura 8. Intervalos para tabular la funcin.
Figura 9. Resultado del programa de Regla Falsa.
1.2. MTODOS ABIERTOS.

En los mtodos cerrados del captulo anterior la raz se encuentra dentro de un intervalo
predeterminado por un lmite inferior y otro superior. La aplicacin repetida de estos mtodos
siempre genera aproximaciones cada vez ms cercanas a la raz. Se dice que tales mtodos
son convergentes porque se acercan progresivamente a la raz a medida que se avanza en el
clculo.
En contraste, los mtodos abiertos descritos en este captulo se basan en frmulas que
requieren nicamente de un solo valor de inicio x o que empiecen con un par de ellos, pero que
no necesariamente encierran la raz. stos, algunas veces divergen o se alejan de la raz
verdadera a medida que se avanza en el clculo. Sin embargo, cuando los mtodos abiertos
convergen, en general lo hacen mucho ms rpido que los mtodos cerrados.
1.2.1.
MTODO DE NEWTON RAPHSON.
Tal vez, de las frmulas para localizar races, la frmula de NewtonRaphson sea la ms
ampliamente utilizada. Si el valor inicial para la raz es xi, entonces se puede trazar una
tangente desde el punto [xi, f (xi)] de la curva. Por lo comn, el punto donde esta tangente
cruza al eje x representa una aproximacin mejorada de la raz. El mtodo de NewtonRaphson
se deduce a partir de esta interpretacin geomtrica (un mtodo alternativo basado en la serie
de Taylor). Se tiene que la primera derivada en x es equivalente a la pendiente:
Que se arregla para obtener:
Figura 10. Frmula de NewtonRaphson.

No hay un criterio general de convergencia para el mtodo de NewtonRaphson. Su
convergencia depende de la naturaleza de la funcin y de la exactitud del valor inicial. La nica
solucin en estos casos es tener un valor inicial que sea suficientemente cercano a la raz.
Para obtener el valor inicial se usan recursos alternativos, tales como las grficas, que
proporcionan mayor claridad en el comportamiento de la solucin. Ante la falta de un criterio
general de convergencia se sugiere el diseo de programas computacionales eficientes que
reconozcan la convergencia lenta o la divergencia.
Aunque en general el mtodo de NewtonRaphson es muy eficiente, hay situaciones donde se

comporta de manera deficiente. En efecto, una pendiente cero [(x) = 0] es un verdadero
desastre, ya que causa una divisin entre cero en la frmula de NewtonRaphson. En forma
grfica, esto significa que la solucin se dispara horizontalmente y jams toca al eje x.
Figura 11. Representacin del mtodo de NewtonRaphson.
1.2.1.1. PROGRAMA Y DIAGRAMA DE FLUJO

Con el fin de hacer prctico y corroborar el resultado arrojado por el mtodo de
biseccin se utilizara los mismos intervalos y la misma funcin.
funcion=4.0 x 330 .0 x2 +70.0 x50.0

#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
float f1(float);
float fp1(float);
int newrap(float, float, float, int, float* );
{
cout<<"Metodo de Newton-Raphson \n"<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de Newton-Raphson \n"<< endl;
cout<<"
<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de Newton-Raphson \n"<< endl;
system("PAUSE");
system("cls");
for(i=1;i<=nr;i++)
{
if (newrap (x1,x2,Eps,nit,&raiz)==1)
{
"<<x1<<endl;
}
else
{
"<<x1<<endl;
}
system("PAUSE");
system("cls");
}
system("PAUSE");
system("cls");
}
int newrap(float xa, float xb, float Eps, int nit, float *res)
{
int con, i;
float y1, y2, yp, x1, x2, x3, Ea;
con=0;
x1=xa;
x2=xb;
x1=(x1+x2)/2.0;
printf(" It
x1
x2 |xi+1 - xi| Ea
y1
y2
yp");
for(i=1; i<=nit;i++)
{
y1=f1(x1);
yp=fp1(x1);
x2=x1-y1/yp;
y2=f1(x2);
if(yp==0.0)
{
cout<<"La derivada es cero "<<endl;
system("PAUSE");
return(0);
}
if(i==1)
printf("\n %2d %8.2f %8.2f \t \t \t %8.2f %8.2f %8.2f", i , x1, x2, y1, y2,
yp);
else
{
Ea=fabs((x2-x3)/x2)*100;
printf("\n %2d %8.2f %8.2f %8.2f %8.2f %8.2f %8.2f %8.2f", i , x1,
x2, fabs(x2-x3), Ea,y1, y2, yp);
}
{
*res=x2;
con=1;
break;
}
else
{
x1=x2;
}
x3=x2;
}
return(con);
}
{
int i;
float x1, dx;
x1 = xa;
dx=(xb-xa)/np;
for (i=1; i<=np ;i++)
{
x[i]=x1;
y[i]= f1(x[i]);
x1=x1+dx;
}
}
{
int i;
for(i=1; i<=n; i++)
{
cout << nomx << "[" <<i<<"]="<<x[i]<<"
"<<nomy<< "["
<<i<<"]="<<y[i]<< endl;
}
}
float f1(float x)
{
}
float fp1(float x)
{
return (12.0*pow(x,2)-60.0*x+70);
}
1.2.1.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN.
Figura 12. Pantalla para tabular la funcin.
Figura 13. Pantalla de resultados.
1.2.2.
MTODO DE LA SECANTE.
Un problema potencial en la implementacin del mtodo de NewtonRaphson es la evaluacin

de la derivada. Aunque esto no es un inconveniente para los polinomios ni para muchas otras
funciones, existen algunas funciones cuyas derivadas en ocasiones resultan muy difciles de
calcular. En dichos casos, la derivada se puede aproximar mediante una diferencia finita
dividida hacia atrs, como en
Esta aproximacin se sustituye en la ecuacin para obtener la siguiente ecuacin iterativa:
Figura 14. La cul es la frmula para el mtodo de la secante.

Observe que el mtodo requiere de dos valores iniciales de x. Sin embargo, debido a que no se
necesita que f(x) cambie de signo entre los valores dados, este mtodo no se clasifica como un
mtodo cerrado.
Figura 15. Representacin grfica del mtodo de la secante.

Esta tcnica es similar a la del mtodo de Newton Raphson en el sentido de que una
aproximacin de la raz se predice extrapolando una tangente de la funcin hasta el eje x. Sin
embargo, el mtodo de la secante usa una diferencia dividida en lugar de una derivada para
estimar la pendiente.

#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
float f1(float);
int sec(float ,float , float , int, float*);
{
cout<<"Metodo de la Secante \n"<<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de la Secante \n"<< endl;
cout<<"
cout<<"xa: Intervalo inicial xb: Intervalo final np: Numero de partes \n" <<endl;
system("PAUSE");
system("cls");
cout<<"
Metodo de la Secante \n"<< endl;
system("PAUSE");
system("cls");
for(i=1;i<=nr;i++)
{
if (sec(x1,x2,Eps,nit,&raiz)==1)
{
Para el intervalo x2: "<<x2<<endl;
}
else
{
Para el intervalo x2: "<<x2<<endl;
system("PAUSE");
system("cls");
}
system("PAUSE");
system("cls");
}
int sec(float xa, float xb, float Eps, int nit, float *res)
{
int i, con;
float x1, x2, x3, x4, Ea, y1, y2, y3, yp;
con=0;
x1=xa;
x2=xb;
printf(" It x1
x2
x3 |xi+1 - xi| Ea
for(i=1;i<=nit;i++)
{
y1=f1(x1);
y1
y2
y3
yp");
y2=f1(x2);
yp=(y2-y1)/(x2-x1);
x3=x2-( y2/yp );
y3=f1(x3);
if(i==1)
{
printf("\n %2d %8.4f %8.4f %8.4f \t \t \t %8.4f %8.4f %8.4f %8.4f", i , x1, x2, x3,
y1, y2, y3, yp);
}
else
{
Ea=fabs((x3-x4)/x3)*100;
printf("\n %2d %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f %8.4f", i ,
x1, x2, x3, fabs(x2-x3), Ea,y1, y2, y3,yp);
{
*res=x3;
con=1;
break;
}
}
x1=x2;
x2=x3;
x4=x3;
}
return(con);
}
{
int i;
float x1, dx;
x1 = xa;
dx=(xb-xa)/np;
for (i=1; i<=np ;i++)
{
x[i]=x1;
y[i]= f1(x[i]);
x1=x1+dx;
}
}
{
int i;
for(i=1; i<=n; i++)
{
cout << nomx << "[" <<i<<"]="<<x[i]<<"

}
"<<nomy<< "[" <<i<<"]="<<y[i]<< endl;
float f1(float x)
{
}
Figura 16. Tabulacion de la funcion.
Figura 17. Resultado del programa de la Secante.
1.3. RACES DE POLINOMIOS

1.3.1.
Mtodo de Lin Bairstow
El mtodo de Bairstow es un mtodo iterativo relacionado de alguna manera con los mtodos
de Mller y de NewtonRaphson. Antes de hacer la descripcin matemtica de ste, recuerde la
forma factorizada de un polinomio.
Con estas consideraciones se puede elaborar un algoritmo para determinar la raz de un
polinomio:
1)
D un valor inicial para la raz x = t.
2)
Divida el polinomio entre el factor x t
3)
Determine si hay un residuo diferente de cero. Si no, el valor inicial es perfecto y la raz
es igual a t. Si existe un residuo, se ajusta el valor inicial en forma sistemtica y se repite el
procedimiento hasta que el residuo desaparezca y se localice la raz. Una vez hecho esto, se
repite el procedimiento totalmente, ahora con el cociente para localizar otra raz.
Por lo general, el mtodo de Bairstow se basa en esta manera de proceder. Por consiguiente,
depende del proceso matemtico de dividir un polinomio entre un factor.
Por ejemplo, el polinomio general:

n (x) = a0 + a1x + a2x2 + + anxn
Se divide entre el factor x t para dar un segundo polinomio que es de un grado menor:
n 1 (x) = b1 + b2x + b3x2 + + bn xn 1
Con un residuo R = b0, donde los coeficientes se calculan por la relacin de recurrencia
bn = an
bi = ai + bi+1t
para i = n 1 a 0
Observe que si t es una raz del polinomio original, el residuo b0 sera igual a cero.
Para permitir la evaluacin de races complejas, el mtodo de Bairstow divide el polinomio entre
un factor cuadrtico x2 rx s. Si esto se hace con la ecuacin, el resultado es un nuevo
polinomio
n 2(x) = b2 + b3x + + bn 1xn3 + bn xn2
Con un residuo
R = b1(x r) + b0
Como con la divisin sinttica normal, se utiliza una relacin de recurrencia simple para realizar
la divisin entre el factor cuadrtico:
bn = an
bn-1 = an-1 + rbn
bi = ai + rbi+1 + sbi+2
para i = n 2 a 0
El factor cuadrtico se introduce para permitir la determinacin de las races complejas. Esto se
relaciona con el hecho de que, si los coeficientes del polinomio original son reales, las races
complejas se presentan en pares conjugados. Si x2 rx s es un divisor exacto del polinomio,
las races complejas pueden determinarse con la frmula cuadrtica. As, el mtodo se reduce
a determinar los valores de r y s que hacen que el factor cuadrtico sea un divisor exacto. En
otras palabras, se buscan los valores que hacen que el residuo sea igual a cero.
La inspeccin de la ecuacin nos lleva a concluir que para que el residuo sea cero, b0 y b1
deben ser cero. Como es improbable que los valores iniciales para evaluar r y s conduzcan a
este resultado, debemos determinar una forma sistemtica para modificar los valores iniciales,
de tal forma que b0 y b1 tiendan a cero. Para lograrlo, el mtodo de Bairstow usa una
estrategia similar a la del mtodo de NewtonRaphson. Como tanto b0 como b1 son funciones
de r y s, se pueden expandir usando una serie de Taylor y se obtiene un sistema de ecuaciones
de derivadas parciales dnde las variables son incrementos de los factores de la aproximacin
del factor cuadrtico, r y
s. Bairstow demostr que las derivadas parciales se obtienen por divisin sinttica de las b en
forma similar a como las b mismas fueron obtenidas:
en = bn
en-1 = bn-1 + ren
ei = bi + rei+1 + sei+2 para i = n 2 a 1
Las derivadas parciales se obtienen por la divisin sinttica de las b y se sustituyen en el
sistema de ecuaciones de derivadas parciales junto con las b para dar el siguiente sistema:
e2r + e3s = b1
e1r + e2s = b0
Estas ecuaciones se resuelven para r y s, las cuales, a su vez, se emplean para mejorar los
valores iniciales de r y s. En cada paso, se estima un error aproximado en r y s:
Cuando ambos errores estimados caen por debajo de un criterio especificado de terminacin
es, los valores de las races se determinan mediante:
En este punto, existen tres posibilidades:

1)
El cociente es un polinomio de tercer grado o mayor. En tal caso, el mtodo de

Bairstow se aplica al cociente para evaluar un nuevo valor de r y s. Los valores
anteriores de r y s pueden servir como valores iniciales en esta aplicacin.
2)
El cociente es cuadrtico. Aqu es posible evaluar directamente las dos races restantes
con la ecuacin.
3)
El cociente es un polinomio de primer grado. En este caso, la raz restante se evala
simplemente como

#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
float f1(float, float, float, float, float);
float linbar(float , float , float, float, float, float);
{
int x;
float a0, a1, a2, a3, Eps, nit;
cout<<"\n Metodo de Lin Baristov"<< endl;
cout<<"\n Inserte la constante de x^3: ";cin>>a3;
cout<<"\n Inserte la constante de x^2: ";cin>>a2;
cout<<"\n Inserte la constante de x: ";cin>>a1;
cout<<"\n Inserte la constante: ";cin>>a0;
f1(x, a0, a1, a2, a3);
cout<<"\n El Polinomio de Tercer Grado es:" <<a3<<"x^3 + "<<a2<<"x^2 +
"<<a1<<"x + "<<a0<<endl;
linbar(Eps, nit, a0, a1, a2, a3);
system("PAUSE");
return EXIT_SUCCESS;
}
float f1(float x, float a0, float a1, float a2, float a3)
{
return ((a3*pow(x,3))+(a2*pow(x,2))+(a1*x)+(a0));
}
float linbar(float Eps, float nit, float a0, float a1, float a2, float a3)
{
int i, n, q, r, s;
float alpha0, alpha1,a[5], b[5], c[5];
a[0]=a0;
a[1]=a1;
a[2]=a2;
a[3]=a3;
alpha0=(a1)/(a2);
alpha1=(a0)/(a2);
if(a2==0)
{
a2=1;
}
n=3;
for(i=0;i<=n;i++)
{
q=i-2;
b[q];
if(q==3)
{
b[3]=0;
}
if(q==2)
{
b[2]=0;
}
b[q]=a[i]-(alpha0*b[i-1])-(alpha1*b[i]);
if(q==1)
{
b[1]=a3-(alpha0*b[2])-(alpha1*b[3]);
}
if(q==0)
{
b[0]=a2-(alpha0*b[1])-(alpha1*b[2]);
}
if(q==-1)
{
b[-1]=a1-(alpha0*b[0])-(alpha1*b[1]);
}
if(q==-2)
{
b[-2]=a0-(alpha0*b[-1])-(alpha1*b[0]);
}
}
for(i=0;i<=n;i++)
{
r=i;
c[r];
}
if(r==2)
{
c[2]=0;
}
if(r==3)
{
c[3]=0;
}
}
Figura 18. El programa pide los coeficientes de la ecuacin de 3 grado.
Figura 19. Resultados del programa.
2. SISTEMAS DE ECUACIONES ALGEBRAICAS LINEALES.
2.1. Eliminacin Gaussiana Simple y Sustitucin hacia atrs

Las ecuaciones algebraicas lineales simultneas que en general se representan como:
Figura 20. Sistema de ecuaciones de nxn.

Donde las a son los coeficientes constantes y las b son los trminos independientes
constantes. La tcnica que se describe en este captulo se conoce como la eliminacin de
Gauss, ya que implica una combinacin de ecuaciones para eliminar las incgnitas. Aunque
ste es uno de los mtodos ms antiguos para resolver ecuaciones lineales simultneas,
contina siendo uno de los algoritmos de mayor importancia, y es la base para resolver
ecuaciones lineales en muchos paquetes de software populares.
La estrategia bsica consiste en multiplicar las ecuaciones por constantes, de tal forma que se
elimine una de las incgnitas cuando se combinen las dos ecuaciones. El resultado es una sola
ecuacin en la que se puede despejar la incgnita restante. Este valor se sustituye en
cualquiera de las ecuaciones originales para calcular la otra variable.
Esta tcnica bsica puede extenderse a sistemas grandes de ecuaciones desarrollando un
esquema sistemtico o algortmico para eliminar incgnitas y sustituir hacia atrs. Aunque tales
tcnicas son muy adecuadas para utilizarlas en computadoras, se requiere de algunas
modificaciones para obtener un algoritmo confiable. En particular, el programa debe evitar la
divisin entre cero. Al mtodo siguiente se le llama eliminacin gaussiana simple, ya que no
evita este problema. En las siguientes secciones se vern algunas caractersticas adicionales
necesarias para obtener un programa de cmputo efectivo.
La primera fase consiste en reducir el conjunto de ecuaciones a un sistema triangular superior.
El paso inicial ser eliminar la primera incgnita, x1, desde la segunda hasta la nsima
ecuacin.
Donde los valores de an,m y bn se calculan como (para k=1,2,3,,N 1 [Tamao del sistema]):
Donde el superndice prima indica que los elementos han cambiado sus valores originales. El
elemento ak,kse denomina el coeficiente o elemento pivote. En particular cuando dicho
elemento es cero hay una divisin entre cero y el mtodo falla, para eso en el Programa de
eliminacin Gaussiana Simple desarrollado en C++ se implementa un algoritmo que
intercambia filas hasta encontrar un elemento diferente de cero, si no se encuentra un valor
diferente a cero el programa indica que el sistema no tiene solucin o hay infinidad de
soluciones. Ya que en dado caso no habra matriz inversa del sistema o su determinante sera
igual a cero, por lo cual se llega a la misma conclusin de que no existe solucin o hay infinidad
de ellas.
Al final de la implementacin del mtodo de Eliminacin Gaussiana Simple, se llega a un
sistema donde se tiene una matriz triangular superior, cmo se muestra a continuacin.
El siguiente paso para resolver el sistema es emplear sustitucin hacia atrs. Este resultado se
puede sustituir hacia atrs en la (n 1) sima ecuacin y despegar xn 1. El procedimiento,
que se repite para evaluar las x restantes, se representa mediante la frmula:
2.1.1. PROGRAMA Y DIAGRAMA DE FLUJO
Para probar el programa se utilizara el siguiente sistema de ecuaciones:
PROGRAMA
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
{
int n, i, j, k;
float aux, a[30][30],x[30];
cout<<"Solucion de Sistema de Ecuaciones Lineales"<<endl;
cout<<"Metodo de Eliminacion de Gauss"<<endl;
cout<<" \n ";
cout<<" \n Numero de Ecuacuaciones: ";
cin>>n;
system("pause");
system("cls");
cout<<" \n ";
cout<<" \n Dame los valores de la matriz: "<<endl;
cout<<" \n "<<endl;
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
cout<< "a[" <<i<<"]["<<j<<"]= ";
cin>> a[i][j];
}
}
system("pause");
system("cls");
for(i=1;i<=n-1;i++)
{
aux=a[i][i];
for(j=i;j<=n+1;j++)
{
a[i][j]=(a[i][j])/aux;
}
for(j=i+1;j<=n;j++)
{
aux=a[j][i];
for(k=i;k<=n+1;k++)
{
a[j][k]=(a[j][k])-((a[i][k])*aux);
}
}
}
x[n]=(a[n][n+1])/(a[n][n]);
for(i=n-1;i>=1;i--)
{
x[i]=a[i][n+1];
k=i+1;
for(j=k;j<=n;j++)
{
x[i]=(x[i])-(a[i][j])*(x[j]);
}
}
cout<<"Solucion"<<endl;
cout<<"\n "<<endl;
for(i=1;i<=n;i++)
{
cout<<"x["<<i<<"]="<<x[i]<<endl;
}
system("PAUSE");
}
2.1.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN.
Figura 21. El programa pregunta el nmero de ecuaciones para generar la matriz.
Figura 22. Pide los valores de la matriz.
Figura 23. Resultados del programa, Siendo acertados los resultados.
2.2. GaussJordn
El mtodo de GaussJordn es una variacin de la eliminacin de Gauss. La principal
diferencia consiste en que cuando una incgnita se elimina en el mtodo de GaussJordn,
sta es eliminada de todas las otras ecuaciones, no slo de las subsecuentes. Adems, todos
los renglones se normalizan al dividirlos entre su elemento pivote. De esta forma, el paso de
eliminacin genera una matriz identidad en vez de una triangular. En consecuencia, no es
necesario usar la sustitucin hacia atrs para obtener la solucin.
Figura 24. Representacin grfica del mtodo de GaussJordn.
Observe la diferencia entre esta tcnica y la de eliminacin de Gauss. El superndice (n)

significa que los elementos del vector del lado derecho se han modificado n veces (en este
caso n = 3).
Aunque la tcnica de GaussJordn y la eliminacin de Gauss podran parecer casi idnticas,
la primera requiere ms trabajo, involucra aproximadamente 50 por ciento ms operaciones
que la eliminacin de Gauss. Por tanto, la eliminacin de Gauss es el mtodo de eliminacin
sencilla que se prefiere para obtener las soluciones de ecuaciones algebraicas lineales.

Se introduce el mismo sistema que se resolvi por Eliminacin Gaussiana Simple para el
programa desarrollado para GaussJordn.
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
{
int n, i, j, k;
float aux, a[30][30],x[30];
cout<<"Metodo de Eliminacion de Gauss-Jordan"<<endl;
cout<<" \n ";
cin>>n;
system("pause");
system("cls");
cout<<" \n ";
cout<<" \n Dame los valores de la matriz: "<<endl;
cout<<" \n "<<endl;
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
cout<< "a[" <<i<<"]["<<j<<"]= ";
cin>> a[i][j];
}
}
system("pause");
system("cls");
for(i=1;i<=n;i++)
{
aux=a[i][i];
for(j=i;j<=n+1;j++)
{
a[i][j]=a[i][j]/aux;
}
for(j=1;j<=n;j++)
{
if(j!=i)
{
aux=a[j][i];
for(k=i;k<=n+1;k++)
{
a[j][k]=(a[j][k])-((a[i][k])*aux);
}
}
}
}
cout<<"Solucion"<<endl;
cout<<"\n "<<endl;
for(i=1;i<=n;i++)
{
x[i]=a[i][n+1];
}
for(i=1;i<=n;i++)
{
x[i];
cout<<"x["<<i<<"]="<<x[i]<<endl;
}
system("PAUSE");
}
Figura 25. El programa pregunta el nmero de ecuaciones para generar la matriz.
Figura 26. Se introducen los datos de la matriz.
Figura 27. Resultado del programa, concuerda con el resultado del mtodo de gauss.
2.3. Matriz Inversa y Multiplicacin de Matrices

Usando el principio la matriz aumentada, se aade en sus columnas a la derecha de una matriz
cuadrada, la matriz identidad, al realizar el mtodo de Gauss Jordn, la matriz del lado del lado
izquierdo resulta ser la identidad y la del lado izquierdo de igual tamao es la inversa. En el
programa se introduce una matriz y se obtiene la inversa, luego se verifica que efectivamente
sta sea su inversa, por lo cual se multiplica por su matriz original y se obtiene la identidad.

PROGRAMA
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
void leematinv(int , float [30][30]);
void escmatinv(int, float [30][30], float [30][30]);
float invmat(int , float [30][30],float [30]);
{
int n;
float a[30][30],b[30][30], x[30];
cout<<"Inversion de Matrices"<<endl;
cin>>n;
cout<<"\n "<<endl;
leematinv(n, a);
system("pause");
system("cls");
invmat(n, a,x);
cout<<"Inversion de Matrices"<<endl;
cout<<"\n "<<endl;
cout<<"La matriz invertida es: "<<endl;
cout<<"\n "<<endl;
escmatinv(n, a , b);
system("PAUSE");
}
float invmat(int n, float a[30][30],float x[30])
{
int i, j, k;
float aux;
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
if(i==j)
{
a[i][n+j]=1.0;
}
else
{
a[i][n+j]=0;
}
}
}
for(i=1;i<=n;i++)
{
aux=a[i][i];
for(j=i;j<=2*n;j++)
{
a[i][j]=a[i][j]/aux;
}
for(j=1;j<=n;j++)
{
if(j!=i)
{
aux=a[j][i];
for(k=i;k<=2*n;k++)
{
a[j][k]=(a[j][k])-((a[i][k])*aux);
}
}
}
}
}
void leematinv(int n, float a[30][30])
{
int i,j;
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
cout<< "a[" <<i<<"]["<<j<<"]= ";
cin>> a[i][j];
}
}
}
void escmatinv(int n, float a[30][30],float b[30][30])
{
int i,j;
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
b[i][j]=a[i][j+n];
cout<< "b[" <<i<<"]["<<j<<"]= "<<b[i][j]<<endl;
}
}
}
Figura 28. Introducir los datos de la matriz que se quiere invertir
Figura 29. Resultados que arroja el programa.
2.4. Mtodo Iterativo de Gauss

Para ste mtodo, se considera el siguiente sistema de ecuaciones.
Y se representa cada una de las variables en trminos de ellas mismas.
Lo cual sugiere el siguiente esquema iterativo de soluciones.
En general se puede escribir como:
El criterio de convergencia que se maneja para ste mtodo es el valor absoluto de la norma
del vector de soluciones de la iteracin anterior menos la norma del vector de soluciones de la
iteracin actual, entre la norma del vector de soluciones de la iteracin anterior.

Para probar el mtodo de Gauss Iterativo implementado en C++, se usa el siguiente sistema
de ecuaciones.
PROGRAMA
#include <cstdlib>
#include <iostream>
#include <math.h>
#include <iomanip>
#include <stdio.h>
{
int i, j, k, con, n, nit;
float Eps, a[30][30],x[30],Xs[30],Xc[30],aux;
cout<<"Metodo de Iterativo de Gauss"<<endl;
system("pause");
system("cls");
cout<<"Inserte el numero de ecuaciones= ";
cin>>n;
cout<<"Inserte el numero iteraciones= ";
cin>>nit;
cout<<"Inserte el Eps= ";
cin>>Eps;
system("pause");
system("cls");
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
cout<<"a["<<i<<"]["<<j<<"]= ";
cin>>a[i][j];
}
}
system("pause");
system("cls");
for(i=1;i<=n;i++)
{
aux=a[i][i];
for(j=1;j<=n;j++)
{
if(i==j)
{
a[i][j]=(a[i][n+1])/aux;
}
else
{
a[i][j]=(-a[i][j])/aux;
}
}
}
for(i=1;i<=n;i++)
{
Xs[i]=a[i][i];
}
for(i=1;i<=nit;i++)
{
for(j=1;j<=n;j++)
{
Xc[j]=a[j][j];
for(k=1;k<=n;k++)
{
if(j!=k)
{
Xc[j]=Xc[j]+(a[j][k]*Xs[k]);
}
}
}
con=0;
for(j=1;j<=n;j++)
{
if(fabs(Xc[j]-Xs[j])>Eps)
{
con=1;
break;
}
}
if(con==0)
{
cout<<"\n Iteraccion: "<<i<<endl;
for(j=1;j<=n;j++)
{
cout<<"\n Converge";
cout<<"\n x["<<j<<"]="<<Xc[j]<<endl;
}
break;
}
else
{
for(j=1;j<=n;j++)
{
Xs[j]=Xc[j];
}
}
}
if(con==1)
{
cout<<"\n No converge";
}
system("PAUSE");
}
Figura 30. Caratula Mtodo Iterativo de Gauss.
Figura 31. Caratula para insertar los valores que ocupa el programa
Figura 32. Se introduce los datos de la matriz del 1 al 3 es el coeficiente a y el cuarto

rengln de la matriz es coeficiente b.
Figura 33. Resultado del programa.
3. SISTEMAS DE ECUACIONES NO LINEALES.

3.1. Mtodo de Newton Raphson para n ecuaciones
Dado el caso general para resolver n ecuaciones no lineales simultneas.
Utilizando la serie de Taylor se puede hacer una aproximacin lineal para una funcin
en un incremento
Si se escribe el sistema en forma matricial, se tiene:
como:
En forma compacta:
Donde J(x) es la matriz de primeras derivada o Jacobiano. Como se quiere encontrar el cero de
la funcin, la aproximacin lineal que se debe resolver es:
Donde las actualizaciones se hacen como:
En el programa desarrollado en C++ la matriz Jacobiana inversa por el vector de funciones, se

define como delta y se resuelve usando Eliminacin Gaussiana Simple.
Se resuelve el siguiente sistema de ecuaciones dado por:
El Jacobiano es:
Y el arreglo de funciones
Considerando como valores iniciales:

El programa se probara con la siguiente funcin, en un intervalo de (1,2)
f ( x , y ) =x3 + x4
PROGRAMA
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <string.h>
#include <conio.h>
#include <time.h>
int main()
{
int i;
float p0, pn, er, tol; // Variables Utilizadas |
printf("\t\tMETODO DE NEWTON-RAPHSON Para encontrar la raiz de la\n");
printf("\t\t funcion F(X)= X^3 + X - 4 , en el intervalo [1,2]\n\n\n");
printf("Ingrese aproximacion inicial P[0] = ");
scanf("%f",&p0);
printf("\nIngrese tolerancia Tol = ");
scanf("%f",&tol);
printf("\n\n\t Pn\t\t\tPorcentaje de Error\n\n");
i=1; // Valores Iniciales
er=1;
while (er>=tol)
{
pn = p0 - (((pow(p0,3)+ p0-4))/(3*(pow(p0,2))-1)); // Clculo de Sucesiones
er=(fabs((pn-p0)/pn)); // Clculo de errores
printf(" P[%d]= %f\t\t",i, pn); // Impresin de resultados
printf("porc. error = %f \n",er);
i=i+1; // Contador n-simo termino
p0=pn; // Redefinicin de Pn
}
printf("\n\nPor lo tanto la RAIZ es = %f", p0); // Impresin Solucin
printf("\n\nEncontrado con %d Iteraciones\n\n", i-1);
system("pause"); //Finaliza y cierra el programa
}
Figura 34. Pantalla del programa en donde se ingresan la tolerancia y la aproximacin

inicial.
4.
AJUSTE DE CURVAS
Es comn que los datos se dan como valores discretos a lo largo de un continuo. Sin embargo,
quizs usted requiera la estimacin de un punto entre valores discretos. Adems, usted puede
necesitar la versin simplificada de una funcin complicada. Una manera de hacerlo es calcular
valores de la funcin en un nmero discreto de valores en el intervalo de inters. Despus, se
obtiene una funcin ms simple para ajustar dichos valores. Estas dos aplicaciones se conocen
como ajuste de curvas.
Existen dos mtodos generales para el ajuste de curvas que se distinguen entre s al
considerar la cantidad de error asociado con los datos. Primero, si los datos exhiben un grado
significativo de error o ruido, la estrategia ser obtener una sola curva que represente la
tendencia general de los datos. Como cualquier dato individual puede ser incorrecto, no se
busca intersecar todos los puntos. En lugar de esto, se construye una curva que siga la
tendencia de los puntos tomados como un grupo. Un procedimiento de este tipo se llama
regresin por mnimos cuadrados.
Segundo, si se sabe que los datos son muy precisos, el procedimiento bsico ser colocar una
curva o una serie de curvas que pasen por cada uno de los puntos en forma directa.
Usualmente tales datos provienen de tablas. Como ejemplos se tienen los valores de la
densidad del agua o la capacidad calorfica de los gases en funcin de la temperatura. La
estimacin de valores entre puntos discretos bien conocidos se llama interpolacin.
4.1. Mnimos Cuadrados

Para entender el procedimiento de mnimos cuadrados para el ajuste de datos con un
polinomio de grado n. Habra que suponer que se ajusta un polinomio de segundo grado o
cuadrtico:
y = a0 + a1x + a2x2 + e
En este caso, la suma de los cuadrados de los residuos es:
Al seguir el procedimiento de la seccin anterior, obtenemos la derivada de la ecuacin con

respecto a cada uno de los coeficientes desconocidos del polinomio
Estas ecuaciones se igualan a cero y se reordenan para desarrollar el siguiente conjunto de

ecuaciones normales:
Donde todas las sumatorias van desde i = 1 hasta n. Observe que las tres ecuaciones
anteriores son lineales y tienen tres incgnitas: a0, a1 y a2. Los coeficientes de las incgnitas
se evalan de manera directa, a partir de los datos observados. En este caso, observamos que
el problema de determinar un polinomio de segundo grado por mnimos cuadrados es
equivalente a resolver un sistema de tres ecuaciones lineales simultneas. En la parte tres se

estudiaron las tcnicas para resolver tales ecuaciones.
El caso bidimensional se extiende con facilidad a un polinomio de nsimo grado como sigue
y = a0 + a1x + a2x2 + + anxn + e
El anlisis anterior se puede extender fcilmente a este caso ms general. As, se reconoce
que la determinacin de los coeficientes de un polinomio de msimo grado es equivalente a
resolver un sistema de m + 1 ecuaciones lineales simultneas. En forma general, para ajustar a
un polinomio de grado n, se tiene:

Para el programa desarrollado en C++ Se ajusta a un polinomio de segundo orden los datos de
la siguiente tabla.
PROGRAMA
#include<iostream>
#include<cstdlib>
#include<stdio.h>
#include<math.h>
#include<iomanip>
void leevec(int, float [30], char [30]);
{
int n, i, j,k,l;
float x[30],y[30], b[30],pol[30][30],aux[30],p[30],v1[30],v2[30],v3[30],a[30],m;
cout<<"Metodo Interpolacion de Lagrange"<<endl<<endl;

cout<<"Ingrese el numero de puntos= ";
cin>>n;
leevec(n,a,"X");
cout<<endl;
leevec(n,y,"Y");
cout<<endl;
for(i=1;i<=n;i++)
{
b[i]=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
b[i]=b[i]*(x[i]-x[j]);
}
}
b[i]=y[i]/b[i];
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
pol[i][j]=0;
}
}
for(i=1;i<=n;i++)
{
k=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
aux[k]=j;
k=k+1;
}
for(j=1;j<=n;j++)
{
v1[j]=0;
v2[j]=0;
}
v1[1]=1;
v1[2]=-m*(aux[1]);
for(j=2;j<=n-1;j++)
{
for(k=1;k<=j;k++)
{
v2[k+1]=v1[k]*(-m*(aux[j]));
}
for(k=1;k<=j;k++)
{
v1[k]=v1[k]+v2[k];
}
for(k=1;k<=j;k++)
{
v2[k]=0;
}
}
}
for(j=1;j<=n;j++)
{
pol[i][j]=v1[j]*b[i];
}
}
for(i=1;i<=n;i++)
{
p[i]=0;
for(j=1;j<=n;j++)
{
p[i]=p[i]+pol[j][i];
}
}
cout<<"El valor del polinomio son:"<<endl;
for(i=1;i<=n;i++)
{
cout<<"p["<<i<<"]="<<p[i]<<endl;
}
system("PAUSE");
}
void leevec(int n, float a[30], char nom[30])
{
int i;
for(i=1;i<=n;i++)
{
cout<<nom<<"["<<i<<"]=";
cin>>a[i];
}
}
Figura 35. Pantalla de inicio para introducir los datos del polinomio y los puntos de
interpolacin.
Figura 36. Resultado del problema con el programa.
4.2. INTERPOLACIN POLINOMIAL DE NEWTON EN DIFERENCIAS

DIVIDIDAS
Como se dijo antes, existe una gran variedad de formas alternativas para expresar una
interpolacin polinomial. El polinomio de interpolacin de Newton en diferencias divididas es
una de las formas ms populares y tiles. Antes de presentar la ecuacin general,
estudiaremos las versiones de primero y segundo grados por su sencilla interpretacin visual.
La forma ms simple de interpolacin consiste en unir dos puntos con una lnea recta. Dicha
tcnica, llamada interpolacin lineal, se ilustra de manera grfica en la figura 4.2. Utilizando
tringulos semejantes.
Reordenndose se tiene:
Que es una frmula de interpolacin lineal. La notacin f1(x) designa que ste es un polinomio
de interpolacin de primer grado. Observe que adems de representar la pendiente de la lnea
que une los puntos, el trmino [f(x1) f(x0)]/(x1 x0) es una aproximacin en diferencia
dividida finita a la primer derivada. En general, cuanto menor sea el intervalo entre los datos,
mejor ser la aproximacin. Esto se debe al hecho de que, conforme el intervalo disminuye,
una funcin continua estar mejor aproximada por una lnea recta.
Figura 37. Ejemplos de interpolacin polinomial: a) de primer grado (lineal) que une dos
puntos, b) de segundo grado (cuadrtica o parablica) que une tres puntos y
c) de tercer grado (cbica) que une cuatro puntos.
Si se tienen tres puntos como datos, stos pueden ajustarse en un polinomio de segundo grado
(tambin conocido como polinomio cuadrtico o parbola). Una forma particularmente
conveniente para ello es
f2(x) = b0 + b1(x x0) + b2(x x0)(x x1)
Observe que aunque la ecuacin parece diferir del polinomio general
f(x) = a0 + a1x + a2x2 + + anxn, las dos ecuaciones son equivalentes.
El anlisis anterior puede generalizarse para ajustar un polinomio de nsimo grado a n+ 1
datos.
El polinomio de nsimo grado es
fn(x) = b0 + b1(x x0) + + bn(x x0) (x x1) (x xn1)
Como se hizo antes con las interpolaciones lineales y cuadrticas, los puntos asociados con
datos se utilizan para evaluar los coeficientes b0, b1,..., bn.
Para un polinomio de nsimo grado se requieren n + 1 puntos: [x0, f(x0)], [x1, f(x1)],..., [xn,
f(xn)]. Se usan estos datos y las siguientes ecuaciones para evaluar los coeficientes:
Donde las evaluaciones de la funcin colocadas entre parntesis son diferencias divididas
finitas. Por ejemplo, la primera diferencia dividida finita en forma general se representa como:
La segunda diferencia dividida finita, que representa la diferencia de las dos primeras
diferencias divididas, se expresa en forma general como:
En forma similar, la nsima diferencia dividida finita es:
Estas diferencias sirven para evaluar los coeficientes en las ecuaciones, los cuales se
sustituirn en la ecuacin para obtener el polinomio de interpolacin
Que se conoce como polinomio de interpolacin de Newton en diferencias divididas.
Figura 38. Representacin grfica de la naturaleza recursiva de las diferencias divididas

finitas.

Para implementar el programa en C++, se consideran los puntos:
PROGRAMA
#include <cstdlib>
#include <iomanip>
#include <iostream>
#include<conio.h>
#include<stdio.h>
#include<math.h>
void leemat(int ,int , float [50][50]);
void leevec (int , float [50], char [5]);
main()
{
int i,j,k,m,n;
float x[50],y[50],a[50],dif[50][50],difdiv[50][50],b[50],pol[50]
[50],p[50],sumx[50],sumxy[20];
cout<<"Diferencias Divididas de Newton";
cout<<"Numero de Puntos= ";cin>>n;
leevec(n,x,"X");
leevec(n,y,"Y");
for(i=1;i<=n;i++)
{
dif[i][1]=y[i];
}
for(i=1;i<=n-1;i++)
{
for(j=1;j<=n-i;j++)
dif[j][i+1]=(dif[j+1][i]-dif[j][i])/(x[i+j]-x[j]);
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
difdiv[i][j]=0.0;
}
for(j=1;j<=n;j++)
{
difdiv[j][1]=y[j];
}
for(i=1;i<=n-1;i++)
{
for(j=1;j<=n-i;j++)
difdiv[j][i+1]=(difdiv[j+1][i]-difdiv[j][i])/(x[i+j]-x[j]);
}
for(i=1;i<=n;i++)
{
b[i]=difdiv[1][i];
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
pol[i][j]=0.0;
}
pol[1][1]=1;
pol[2][1]=1;
pol[2][2]=-x[1];
for(i=2;i<=n;i++)
{
for(j=1;j<=i;j++)
{
pol[i+1][j+1]=pol[i][j]*(-x[i]);
}
for(j=1;j<=i;j++)
{
pol[i+1][j]=pol[i+1][j]+pol[i][j];
}
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
pol[i][j]=pol[i][j]*b[i];
}
}
for(i=1;i<=n-1;i++)
{
k=2;
p[n+1-i]=pol[i][1];
for(j=i+1;j<=n;j++)
{
p[n+1-i]=p[n+1-i]+pol[j][k];
k=k+1;
}
}
p[1]=pol[n][1];
for(i=1;i<=n;i++)
{
cout<<"\n";
cout<<"["<<i<<"]="<<x[i];
}
}
void leevec(int n, float X[20], char nom[5])
{
int i;
for(i=1;i<=n;i++)
{
cout<<nom<<"["<<i<<"]=";
cin>>X[i];
}
}
Figura 39. Pantalla de inicio para ingresar los valores de los puntos.
Figura 40. Pantalla que arroja los resultados del programa.
4.3. Polinomios de Interpolacin de LaGrange
El polinomio de interpolacin de Lagrange es simplemente una reformulacin del polinomio de

Newton que evita el clculo de las diferencias divididas, y se representa de manera concisa
como
Donde
Donde designa el producto de. Por ejemplo, la versin lineal (n = 1)
Y la versin de segundo grado es
El razonamiento detrs de la formulacin de LaGrange se comprende directamente al darse

cuenta de que cada trmino Li(x) ser 1 en x = xi y 0 en todos los otros puntos. De esta forma,
cada producto Li(x) f (xi) toma el valor de f (xi) en el punto xi. En consecuencia, la sumatoria de
todos los productos en la ecuacin (61) es el nico polinomio de nsimo grado que pasa
exactamente a travs de todos los n + 1 puntos, que se tienen como datos.

Para el programa en C++. Se comprueban los mismos puntos que se emplearon para encontrar
el polinomio de diferencias divididas de Newton, para el polinomio de interpolacin de
Lagrange. ste tiene que ser exactamente el mismo.
PROGRAMA
#include<iostream>
#include<cstdlib>
#include<stdio.h>
#include<math.h>
#include<iomanip>
float interlagr(int ,float [30],float [30], float [30],float [30]);
{
int n, i, j,k;
float r[30],y[30], b[30],p[30],x[30],pol[30]
[30],aux[30],v1[30],v2[30],v3[30];
cout<<"Metodo Interpolacion de Lagrange"<<endl<<endl;
cout<<"Ingrese el numero de puntos= ";
cin>>n;
leevec(n,r,"X");
cout<<endl;
leevec(n,y,"Y");
cout<<endl;
interlagr(n, x, y, b,p);
for(i=1;i<=n;i++)
{
b[i]=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
b[i]=b[i]*(x[i]-x[j]);
}
}
b[i]=y[i]/b[i];
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
pol[i][j]=0;
}
}
for(i=1;i<=n;i++)
{
k=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
aux[k]=j;
k=k+1;
}
for(j=1;j<=n;j++)
{
v1[j]=0;
v2[j]=0;
}
v1[1]=1;
v1[2]=-x[aux[1]];
for(j=2;j<=n-1;j++)
{
for(k=1;k<=j;k++)
{
v2[k+1]=v1[k]*(-
x[aux[j]]);
for(k=1;k<=j;k++)
{
v1[k]=v1[k]+v2[k];
}
for(k=1;k<=j;k++)
{
}
v2[k]=0;
}
for(j=1;j<=n;j++)
{
}
}
for(i=1;i<=n;i++)
{
p[i]=0;
for(j=1;j<=n;j++)
{
}
}
for(i=1;i<=n;i++)
{
cout<<"p["<<i<<"]="<<p[i]<<endl;
}
system("PAUSE");
}
void leevec(int n, float r[30], char nom[30])

{
int i;
for(i=1;i<=n;i++)
{
cout<<nom<<"["<<i<<"]=";
cin>>r[i];
}
}
float interlagr(int n, int m,float x[30],float y[30], float

b[30],float p[30])
{
int
i,j,k;
float pol[30][30],aux[30],v1[30],v2[30],v3[30];
for(i=1;i<=n;i++)
{
b[i]=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
b[i]=b[i]*(x[i]-x[j]);
}
b[i]=y[i]/b[i];
}
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
pol[i][j]=0;
}
}
for(i=1;i<=n;i++)
{
k=1;
for(j=1;j<=n;j++)
{
if(i!=j)
{
aux[k]=j;
k=k+1;
}
for(j=1;j<=n;j++)
{
v1[j]=0;
v2[j]=0;
}
v1[1]=1;
v1[2]=-x[aux[1]];
for(j=2;j<=n-1;j++)
{
for(k=1;k<=j;k++)
{
v2[k+1]=v1[k]*(-
x[aux[j]]);
}
for(k=1;k<=j;k++)
{
v1[k]=v1[k]+v2[k];
}
for(k=1;k<=j;k++)
{
}
v2[k]=0;
}
for(j=1;j<=n;j++)
{
}
}
for(i=1;i<=n;i++)
{
p[i]=0;
for(j=1;j<=n;j++)
{
}
}
for(i=1;i<=n;i++)
{
cout<<"p["<<i<<"]="<<p[i]<<endl;
}
Figura 41. PANTALLA PARA INTRODUCCIR LOS VALORES DE LA FUNCION
Figura 42. Caratula de resultados.
Metodo de Interpolacion de Lagrange

Numero de datos de tu lista 6 Introduce los valores de x[i] x[0]= 0
x[1]= 1
x[2]= 2
x[3]= 3
x[4]= 4
x[5]= 5
Introduce los valores de f(x)
fx(0)= 2.10 fx(1)= 7.70 fx(2)= 13.60 fx(3)= 27.20 fx(4)= 40.90 fx(5)= 61.10
/n El polinomio de Lagrange :
(0.241666x^5) + (3.04167x^4) + (13.4417x^3) + (22.5083x^2) + (17.4667x^1) + (2.1x^0) +
Process exited after 29.03 seconds with return value 0 Presione una tecla para continuar . . .
Figura 4.5 Programa de Interpolacin de Lagrange para un polinomio de grado n1 puntos

dados
V.
ECUACIONES DIFERENCIALES ORDINARIAS.
Sea una ecuacin diferencial ordinaria de la forma
(65)
Y sea
(66)
La ecuacin de pendiente ordenada al origen. De acuerdo con esta ecuacin, la pendiente

estimada f se usa para extrapolar desde un valor anterior yi a un nuevo valor yi+1 en una
distancia h (figura 5.1).
Nuevo valor = valor anterior + pendiente tamao de paso
Figura 5.1 Ilustracin grfica del mtodo de un paso.
Esta frmula se aplica paso a paso para calcular un valor posterior y, por lo tanto, para trazar la
trayectoria de la solucin. Todos los mtodos de un paso que se expresen de esta forma
general, tan slo van a diferir en la manera en la que se estima la pendiente. En otras palabras,
se toma la pendiente al inicio del intervalo como una aproximacin de la pendiente promedio
sobre todo el intervalo. Tal procedimiento, llamado mtodo de Euler. Despus se revisa otro
mtodo de un paso que emplean otras formas de estimar la pendiente que dan como resultado
predicciones ms exactas. Todas estas tcnicas en general se conocen como mtodos de
RungeKutta.
Mtodo de Euler
La primera derivada ofrece una estimacin directa de la pendiente en xi
c/ = (xi, yi) (67)
Donde (xi, yi) es la ecuacin diferencial evaluada en xi y yi. La estimacin se sustituye en la

ecuacin (25.1):
y i+1 = yi + (xi, yi)h (68)

Esta frmula se conoce como mtodo de Euler (o de EulerCauchy o de punto pendiente). Se
predice un nuevo valor de y usando la pendiente (igual a la primera derivada en el valor original
de
x) para extrapolar linealmente sobre el tamao de paso h.
Para el programa implementado en C++ la ecuacin diferencial ordinaria de primer orden se

escribe en la primera librera #define como.
fp(x,y)(exp(x)*2*cos(2*x)exp(x)*sin(2*x))
Como el programa proporciona tambin los puntos de la solucin original, se define tambin
f(x,y)(sin(2*x)*exp(x))
En las directivas #define, la cul es la solucin a dicha EDO. El programa muestra en forma
tabulada los resultados paso a paso.
Programa Mtodo de Euler

*Nota: Introducir datos en formato decimal (sin indicar operaciones. Ejemplo: 6.2832, en vez de
2*PI o 2*3.1416)
Limite inferior del intervalo (ti):0 Limite superior del intervalo (tf):6.2832 Condicin inicial (xi):0
Condicin inicial (yi):0 Ancho de intervalo (h):0.3
xi
yi
ytrue
f(x,y)
+0.000000+0.000000+0.000000+2.000000
+0.300000+0.600000+0.418297+0.804550
+0.600000+0.841365+0.5115140.113782
+0.900000+0.807230+0.3959370.580684
+1.200000+0.633025+0.2034460.647643
+1.500000+0.438732+0.0314880.473282
+1.800000+0.2967480.0731480.223318
+2.100000+0.2297520.1067300.013341
+2.400000+0.2257500.090370+0.106245
+2.700000+0.2576240.051934+0.137244
+3.000000+0.2987970.013911+0.109519
+3.300000+0.331652+0.011491+0.058605
+3.600000+0.349234+0.021686+0.011559
+3.900000+0.352701+0.0202120.018028
+4.200000+0.347293+0.0128150.028389
+4.500000+0.338776+0.0045780.024822
+4.800000+0.3313300.0014350.014773
+5.100000+0.3268980.0042670.004442
+5.400000+0.3255650.004430+0.002675
+5.700000+0.3263680.003076+0.005709
+6.000000+0.3280800.001330+0.005513
+6.300000+0.329735+0.000062+0.003609
Figura 5.2 Programa de Euler.
Mtodos de RungeKutta de cuarto orden

El ms popular de los mtodos RK es el de cuarto orden. Como en el caso de los
procedimientos de segundo orden, hay un nmero infinito de versiones. La siguiente, es la
forma comnmente usada y, por lo tanto, le llamamos mtodo clsico RK de cuarto orden:
(69)
Donde
(70)
Observe que con las EDO que estn en funcin slo de x, el mtodo RK clsico de cuarto
orden es similar a la regla de Simpson 1/3. Adems, el mtodo RK de cuarto orden tiene
similitud con el procedimiento de Heun en cuanto a que se usan mltiples estimaciones de la
pendiente para obtener una mejor pendiente promedio en el intervalo. Como se muestra en la
figura 5.2, cada una de las k representa una pendiente.
Figura 5.3. Representacin grfica de las pendientes estimadas empleadas en el mtodo RK

de cuarto Orden
Para RK de 4to orden, se usa la misma funcin que se us para Euler, con el mismo paso y
condiciones iniciales.
Programa Mtodo de RungeKutta de 4to orden

*Nota: Introducir datos en formato decimal (sin indicar operaciones. Ejemplo: 6.2832, en vez de
2*PI o 2*3.1416)
Limite inferior del intervalo (ti):0 Limite superior del intervalo (tf):6.2832 Condicin inicial (xi):0
Condicin inicial (yi):0 Ancho de intervalo (h):0.3
xi
yi
y true k1
k2
k3
k4
Pend tot
+0.000000+0.000000+0.000000+2.000000+1.390175+1.390175+0.804550+0.418262
+0.300000+0.418262+0.418297+0.804550+0.293241+0.2932410.113782+0.093187
+0.600000+0.511449+0.5115140.1137820.4043550.4043550.5806840.115594
+0.900000+0.395855+0.3959370.5806840.6553990.6553990.6476430.192496
+1.200000+0.203359+0.2034460.6476430.5795380.5795380.4732820.171954
+1.500000+0.031405+0.0314880.4732820.3489960.3489960.2233180.104629
+1.8000000.0732250.0731480.2233180.1087110.1087110.0133410.033575
+2.1000000.1068000.1067300.013341+0.058595+0.058595+0.106245+0.016364
+2.4000000.0904350.090370+0.106245+0.131315+0.131315+0.137244+0.038438
+2.7000000.0519980.051934+0.137244+0.128421+0.128421+0.109519+0.038022
+3.0000000.0139760.013911+0.109519+0.084972+0.084972+0.058605+0.025401
+3.300000+0.011425+0.011491+0.058605+0.033428+0.033428+0.011559+0.010194
+3.600000+0.021619+0.021686+0.0115590.0057550.0057550.0180280.001475
+3.900000+0.020144+0.0202120.0180280.0253840.0253840.0283890.007398
+4.200000+0.012747+0.0128150.0283890.0278820.0278820.0248220.008237
+4.500000+0.004510+0.0045780.0248220.0201650.0201650.0147730.006013
+4.8000000.0015030.0014350.0147730.0093560.0093560.0044420.002832
+5.1000000.0043350.0042670.0044420.0003750.000375+0.0026750.000163
+5.4000000.0044980.004430+0.002675+0.004677+0.004677+0.005709+0.001355
+5.7000000.0031440.003076+0.005709+0.005924+0.005924+0.005513+0.001746
+6.0000000.0013980.001330+0.005513+0.004678+0.004678+0.003609+0.001392
+6.3000000.000006+0.000062+0.003609+0.002469+0.002469+0.001387+0.000744
Figura 5.4 Programa de Runge Kutta de 4to. Orden.
Comparacin Grfica de los Mtodos de Runge Kutta en MATLAB hasta 5to.

Orden.
Se modific el cdigo original para aplicar el mtodo por separado, para ejecutar un script
(exp_sin_b_2.m) por separado llamando a las funciones originales, en ste caso el mtodo de
Euler, cmo su funcin es la primera en ser llamada, es en el que se grafica la funcin original,
ese cdigo se omite en las funciones de los dems mtodos. Los archivos .m por separado,
estn programados para proporcionar la grfica del mtodo en comparacin con la funcin
original, que se introduzca su derivada, cada funcin recibe dicha funcin y calcula su primitiva
para ser comparada; adems se proporcionan los valores en forma tabulada.
Nota: En los archivos .m originales, hay vario cdigo en comentarios, la mayora es significativo
y sin gran relevancia, ya que se fue perfeccionando poco a poco hasta corregir errores de
programacin y se omite en ste documento. Hay otro cdigo que se puede usar para
proporcionar valores que se omiten en las tabulaciones. Al final de ste texto, se agregan los
cdigos de las funciones utilizadas.
Cdigos del Script para la funcin:

f(x) = e-x sin b x
f(x) = e-x (b cos(b x) - sin(b x))
Para b=2
syms x f;
f='exp(-x)*(2*cos(2*x)-sin(2*x))'; ti=0;
tf=2*pi; Xi=0; Yi=0; h=0.3;
Euler_Method(f,ti,tf,Xi,Yi,h); Heun_Method(f,ti,tf,Xi,Yi,h); Mid_Method(f,ti,tf,Xi,Yi,h);
Runge_Kutta2_Ralston(f,ti,tf,Xi,Yi,h); Runge_Kutta3(f,ti,tf,Xi,Yi,h); Runge_Kutta4(f,ti,tf,Xi,Yi,h);
Runge_Kutta5_Butcher(f,ti,tf,Xi,Yi,h);
title('Comparacin de los mtodos de Runge Kutta con la funcin original');
xlabel('Tiempo t (seg)');
ylabel('f(t)');
legend('Euler','Solucin verdadera','Heun','Punto Medio','Ralston','Runge Kutta 3', 'Runge Kutta
4','Runge Kutta 5 (Butcher)');
1.2
0.8
0.6
Comparacin de los mtodos de Runge Kutta con la funcin original
Euler
Solucin verdadera Heun
Punto Medio Ralston
Runge Kutta 3
Runge Kutta 4
Runge Kutta 5 (Butcher)
0.4
0.2
-0.2
Tiempo t (seg)
Figura 5.5 Comparacin grafica de los mtodos de Runge Kutta, para la EDO.
f(x) = e-x (2 cos(b x) - sin(2 x))
Salida en la ventana de comandos
primitiva = sin(2*x)*exp(-x)
i
Y_Verdadero Y_Euler
Pendiente
Y+1_Euler
B=
1.0000 0.3000 0.4183 0.6000 2.0000 0.6000

2.0000 0.6000 0.5115 0.8414 0.8045 0.8414
3.0000 0.9000 0.3959 0.8072 -0.1138
0.8072
4.0000 1.2000 0.2034 0.6330 -0.5807
0.6330
5.0000 1.5000 0.0315 0.4387 -0.6476
0.4387
6.0000 1.8000 -0.0731
0.2967 -0.4733
0.2967
7.0000 2.1000 -0.1067
0.2298 -0.2233
0.2298
8.0000 2.4000 -0.0904
0.2257 -0.0133
0.2257
9.0000 2.7000 -0.0519
0.2576 0.1062 0.2576
10.0000
3.0000 -0.0139
0.2988 0.1372 0.2988
11.0000
3.3000 0.0115 0.3317 0.1095 0.3317
12.0000
3.6000 0.0217 0.3492 0.0586 0.3492
13.0000
3.9000 0.0202 0.3527 0.0116 0.3527
14.0000
4.2000 0.0128 0.3473 -0.0180
0.3473
15.0000
4.5000 0.0046 0.3388 -0.0284
0.3388
16.0000
4.8000 -0.0014
0.3313 -0.0248
0.3313
17.0000
5.1000 -0.0043
0.3269 -0.0148
0.3269
18.0000
5.4000 -0.0044
0.3256 -0.0044
0.3256
19.0000
5.7000 -0.0031
0.3264 0.0027 0.3264
20.0000
6.0000 -0.0013
0.3281 0.0057 0.3281
x
B
Y_Verdadero Y_Heun
=
0
0.3000 0.4183 0.4207

0.6000 0.5115 0.5243
0.9000 0.3959 0.4201
1.2000 0.2034 0.2359
1.5000 0.0315 0.0677
1.8000 -0.0731
-0.0368
2.1000 -0.1067
-0.0722
2.4000 -0.0904
-0.0583
2.7000 -0.0519
-0.0218
3.0000 -0.0139
0.0152
3.3000 0.0115 0.0404

3.6000 0.0217 0.0510
3.9000 0.0202 0.0500
4.2000 0.0128 0.0430
4.5000 0.0046 0.0351

4.8000 -0.0014
0.0291
5.1000 -0.0043
0.0262
5.4000 -0.0044
0.0260
5.7000 -0.0031
0.0272
6.0000 -0.0013
0.0289
B=
Y_Verdadero
Y_Punto_Medio
0
0.3000 0.4183 0.4171
0.6000 0.5115 0.5050
0.9000 0.3959 0.3837
1.2000 0.2034 0.1871
1.5000 0.0315 0.0132
1.8000 -0.0731
-0.0915
2.1000 -0.1067
-0.1241
2.4000 -0.0904
-0.1065
2.7000 -0.0519
-0.0671
3.0000 -0.0139
-0.0286
3.3000 0.0115 -0.0031

3.6000 0.0217 0.0069
3.9000 0.0202 0.0052
4.2000 0.0128 -0.0024
4.5000 0.0046 -0.0108

4.8000 -0.0014
-0.0168
5.1000 -0.0043
-0.0196
5.4000 -0.0044
-0.0197
5.7000 -0.0031
-0.0183
6.0000 -0.0013
-0.0166
x
B
Y_Verdadero Y_Runge_Kutta2_Ralston
=
0
0.3000 0.4183 0.4181

0.6000 0.5115 0.5137
0.9000 0.3959 0.4011
1.2000 0.2034 0.2109
1.5000 0.0315 0.0402
1.8000 -0.0731
-0.0642
2.1000 -0.1067
-0.0982
2.4000 -0.0904
-0.0825
2.7000 -0.0519
-0.0446
3.0000 -0.0139
-0.0068
3.3000 0.0115 0.0185

3.6000 0.0217 0.0287
3.9000 0.0202 0.0274
4.2000 0.0128 0.0201
4.5000 0.0046 0.0119
4.8000 -0.0014
0.0059
5.1000 -0.0043
0.0031
5.4000 -0.0044
0.0029
5.7000 -0.0031
0.0042
6.0000 -0.0013
0.0060
B=
Y_Verdadero
Y_Runge_Kutta3
0
0.3000 0.4183 0.4183
0.6000 0.5115 0.5114
0.9000 0.3959 0.3959
1.2000 0.2034 0.2034
1.5000 0.0315 0.0314
1.8000 -0.0731
-0.0732
2.1000 -0.1067
-0.1068
2.4000 -0.0904
-0.0904
2.7000 -0.0519
-0.0520
3.0000 -0.0139
-0.0140
3.3000 0.0115 0.0114

3.6000 0.0217 0.0216
3.9000 0.0202 0.0201

4.2000 0.0128 0.0127
4.5000 0.0046 0.0045
4.8000 -0.0014
-0.0015
5.1000 -0.0043
-0.0043
5.4000 -0.0044
-0.0045
5.7000 -0.0031
-0.0031
6.0000 -0.0013
B=
Y_Verdadero
Y_Runge_Kutta4
0
0.3000 0.4183 0.4183
0.6000 0.5115 0.5114
0.9000 0.3959 0.3959
1.2000 0.2034 0.2034
-0.0014
1.5000 0.0315 0.0314

1.8000 -0.0731
-0.0732
2.1000 -0.1067
-0.1068
2.4000 -0.0904
-0.0904
2.7000 -0.0519
-0.0520
3.0000 -0.0139
-0.0140
3.3000 0.0115 0.0114

3.6000 0.0217 0.0216
3.9000 0.0202 0.0201
4.2000 0.0128 0.0127
4.5000 0.0046 0.0045
4.8000 -0.0014
-0.0015
5.1000 -0.0043
-0.0043
5.4000 -0.0044
-0.0045
5.7000 -0.0031
6.0000 -0.0013
-0.0014
B=
Y_Verdadero
Y_Runge_Kutta5_Butcher
-0.0031
0
0.3000 0.4183 0.4183
0.6000 0.5115 0.5115
0.9000 0.3959 0.3959
1.2000 0.2034 0.2034
1.5000 0.0315 0.0315
1.8000 -0.0731
-0.0731
2.1000 -0.1067
-0.1067
2.4000 -0.0904
-0.0904
2.7000 -0.0519
-0.0519
3.0000 -0.0139
-0.0139
3.3000 0.0115 0.0115

3.6000 0.0217 0.0217
3.9000 0.0202 0.0202
4.2000 0.0128 0.0128
4.5000 0.0046 0.0046
4.8000 -0.0014
-0.0014
5.1000 -0.0043
-0.0043
5.4000 -0.0044
-0.0044
5.7000 -0.0031
-0.0031
6.0000 -0.0013
-0.0013
Cdigo de las funciones en MATLAB de los mtodos usados
Cdigo para el Mtodo de Euler:
function Euler_Method(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

x1=zeros(n,1); y1=zeros(n,1); x1=[ti:h:tf]; y1(1)=Yi;
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y=y1(k); fp(k+1)=eval(f);
y2(k+1)=y1(k)+(fp(k+1))*h; y1(k+1)=y2(k+1);
A(k+1,1)=k;
end
plot(x1,y1,'-r') primitiva=int(sym(f)) Q=length(x1);
for u=1:Q
end
x=x1(u); y3(u)=eval(primitiva);
hold on plot(x1,y3,'-b')
disp(' i
Y_Verdadero Y_Euler
B=[A, x1', y3', y1, fp', y2']

end
Cdigo para el Mtodo de Heun:
function Heun_Method(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

Pendiente Y+1_Euler')
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y2(k+1)=y1(k)+(fp(k+1))*h; y=y2(k+1);
x=x1(k+1); fp2(k+1)=eval(f);
feval(k+1)=(fp(k+1)+fp2(k+1))/2; y1(k+1)=y1(k)+feval(k+1)*h; A(k+1,1)=k;
end
plot(x1,y1,'-g') primitiva=int(sym(f)); Q=length(x1);
for u=1:Q
end
disp(' x
Y_Verdadero Y_Heun') B=[x1', y4', y1]
end
Cdigo para el Mtodo de Punto Medio:
function Mid_Method(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

x1=zeros(n,1); y1=zeros(n,1); x1=[ti:h:tf];
y1(1)=Yi;
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y2(k+1)=y1(k)+((fp(k+1))*h)/2; y=y2(k+1);
x=x1(k)+h/2; fp2(k+1)=eval(f); y1(k+1)=y1(k)+(fp2(k+1))*h; A(k+1,1)=k;
end
plot(x1,y1,'-m') primitiva=int(sym(f)); Q=length(x1);
for u=1:Q
end
disp(' x
Y_Verdadero Y_Punto_Medio') B=[x1', y4', y1]
end
Cdigo para el Mtodo de Runge Kutta de 2do. Orden (Ralston):
function Runge_Kutta2_Ralston(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y=y1(k); k1(k+1)=eval(f);
y2(k+1)=y1(k)+((k1(k+1))*3*h)/4; y=y2(k+1);
x=x1(k)+3*h/4; k2(k+1)=eval(f);
y1(k+1)=y1(k)+((1/3)*k1(k+1)+(2/3)*(k2(k+1)))*h; A(k+1,1)=k;
end
plot(x1,y1,'-c') primitiva=int(sym(f)); Q=length(x1);
for u=1:Q
end
disp(' x
Y_Verdadero Y_Runge_Kutta2_Ralston') B=[x1', y4', y1]
end
Cdigo para el Mtodo de Runge Kutta de 3er. Orden:
function Runge_Kutta3(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

x1=zeros(n,1); y1=zeros(n,1); x1=[ti:h:tf]
y1(1)=Yi;
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y2(k+1)=y1(k)+((k1(k+1))*h)/2; y=y2(k+1);
x=x1(k)+h/2; k2(k+1)=eval(f);
y3(k+1)=y1(k)-h*k1(k+1)+2*h*k2(k+1); y=y3(k+1);
x=x1(k)+h; k3(k+1)=eval(f);
y1(k+1)=y1(k)+(1/6)*(k1(k+1)+4*(k2(k+1))+(k3(k+1)))*h; A(k+1,1)=k;
end plot(x1,y1,'ok')
primitiva=int(sym(f)); Q=length(x1);
for u=1:Q
end
disp(' x
Y_Verdadero Y_Runge_Kutta3') B=[x1', y4', y1]
end
Cdigo para el Mtodo de Runge Kutta de 4to. Orden:
function Runge_Kutta4(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

y1(1)=Yi;
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y2(k+1)=y1(k)+((k1(k+1))*h)/2; y=y2(k+1);
end
x=x1(k)+h/2; k2(k+1)=eval(f); y3(k+1)=y1(k)+((k2(k+1))*h)/2; y=y3(k+1);

x=x1(k)+h/2; k3(k+1)=eval(f); y4(k+1)=y1(k)+(k3(k+1))*h; y=y4(k+1);
x=x1(k)+h; k4(k+1)=eval(f);
y1(k+1)=y1(k)+(1/6)*(k1(k+1)+2*(k2(k+1))+2*(k3(k+1))+k4(k+1))*h; A(k+1,1)=k;
plot(x1,y1,'*y') primitiva=int(sym(f)) Q=length(x1);

for u=1:Q
end
disp(' x
Y_Verdadero Y_Runge_Kutta4') B=[x1', y4', y1]
end
Cdigo para el Mtodo de Runge Kutta de 5to. Orden (Butcher):
function Runge_Kutta5_Butcher(f,ti,tf,Xi,Yi,h) n=(tf-ti)/h;

y1(1)=Yi;
y2(1)=Yi;
A(1)=0;
for k=1:n
x=x1(k);
y2(k+1)=y1(k)+((k1(k+1))*h)/4; y=y2(k+1);
x=x1(k)+h/4; k2(k+1)=eval(f);
y3(k+1)=y1(k)+((k1(k+1))*h)/8+((k2(k+1))*h)/8; y=y3(k+1);
x=x1(k)+h/4; k3(k+1)=eval(f);
y4(k+1)=y1(k)-((k2(k+1))*h)/2+(k3(k+1))*h; y=y4(k+1);
x=x1(k)+h/2; k4(k+1)=eval(f);
y5(k+1)=y1(k)+(3/16)*(k1(k+1))*h+(9/16)*h*(k4(k+1)); y=y5(k+1);
x=x1(k)+(3/4)*h;
k5(k+1)=eval(f);
y6(k+1)=y1(k)-(3/7)*(k1(k+1))*h+(2/7)*(k2(k+1))*h+(12/7)*(k3(k+1))*h(12/7)*(k4(k+1))*h+(8/7)*(k5(k+1))*h;
y=y6(k+1); x=x1(k)+h; k6(k+1)=eval(f);
y1(k+1)=y1(k)+(1/90)*(7*(k1(k+1))+32*(k3(k+1))+12*(k4(k+1))+32*(k5(k+1))+
7*(k6(k+1)))*h; A(k+1,1)=k;
end plot(x1,y1,'xk')
primitiva=int(sym(f)) Q=length(x1);
for u=1:Q
end
disp(' x
end
Y_Verdadero Y_Runge_Kutta5_Butcher') B=[x1', y4', y1]

Trabajo de Metodos

Cargado por

Copyright:

Formatos disponibles

Trabajo de Metodos

Cargado por

Información del documento

Derechos de autor

Formatos disponibles

Compartir este documento

Compartir o incrustar documentos

Opciones para compartir

¿Le pareció útil este documento?

¿Este contenido es inapropiado?

Copyright:

Formatos disponibles

Trabajo de Metodos

Cargado por

Copyright:

Formatos disponibles

INSTITUTO POLITECNICO NACIONAL

Escuela Superior de Ingeniera Mecnica y

[PROGRAMACIN Y MTODOS NUMRICOS]

ALUMNO: CARLOS ALFREDO HARRISON ROJAS

PROGRAMA DE TABULACION Y DE BISECCION DE UN POLINOMIO........4

CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN....................7

METODO DE REGLA FALSA........................................................................11

PROGRAMA DE TABULACIN Y DEL MTODO DE REGLA FALSA..........12

CAPTURA DE PANTALLA DEL PROGRAMA EN EJECUCIN....................15

PROGRAMA Y DIAGRAMA DE FLUJO.....................................................18

CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN..................21

PROGRAMA Y DIAGRAMA DE FLUJO.....................................................23

CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN..................26

PROGRAMA Y DIAGRAMA DE FLUJO.....................................................30

CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN..................32

SISTEMAS DE ECUACIONES ALGEBRAICAS LINEALES.........................................33

Eliminacin Gaussiana Simple y Sustitucin hacia atrs..............................33

Figura 1. Mtodo de biseccin

Donde xr nuevo es la raz en la iteracin actual y xr anterior es el valor de la raz en la iteracin

1.1.1.1. PROGRAMA DE TABULACION Y DE BISECCION DE UN POLINOMIO

funcion=4.0 x 330 .0 x2 +70.0 x50.0

1.1.1.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN

Figura 2. El programa solicita intervalos y en cuantas partes hace la tabulacin.

Figura 3. Ventana de resultados de la tabulacin del polinomio.

Figura 4. Ya con la funcin tabulada podemos encontrar un intervalo menor para

Figura 5. Resultados del programa el mtodo converge en 16 iteraciones y nos

METODO DE REGLA FALSA

1.1.2.1. PROGRAMA DE TABULACIN Y DEL MTODO DE REGLA FALSA

funcion=4.0 x 330 .0 x2 +70.0 x50.0

1.1.2.2. CAPTURA DE PANTALLA DEL PROGRAMA EN EJECUCIN

Figura 8. Intervalos para tabular la funcin.

Figura 9. Resultado del programa de Regla Falsa.

1.2. MTODOS ABIERTOS.

MTODO DE NEWTON RAPHSON.

Que se arregla para obtener:

Figura 10. Frmula de NewtonRaphson.

Aunque en general el mtodo de NewtonRaphson es muy eficiente, hay situaciones donde se

Figura 11. Representacin del mtodo de NewtonRaphson.

1.2.1.1. PROGRAMA Y DIAGRAMA DE FLUJO

funcion=4.0 x 330 .0 x2 +70.0 x50.0

1.2.1.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN.

Figura 12. Pantalla para tabular la funcin.

Figura 13. Pantalla de resultados.

Un problema potencial en la implementacin del mtodo de NewtonRaphson es la evaluacin

Esta aproximacin se sustituye en la ecuacin para obtener la siguiente ecuacin iterativa:

Figura 14. La cul es la frmula para el mtodo de la secante.

Figura 15. Representacin grfica del mtodo de la secante.

1.2.2.1. PROGRAMA Y DIAGRAMA DE FLUJO

cout << nomx << "[" <<i<<"]="<<x[i]<<"

"<<nomy<< "[" <<i<<"]="<<y[i]<< endl;

1.2.2.2. CAPTURAS DE PANTALLA DEL PROGRAMA EN EJECUCIN.

Figura 16. Tabulacion de la funcion.

Figura 17. Resultado del programa de la Secante.

1.3. RACES DE POLINOMIOS

D un valor inicial para la raz x = t.

Divida el polinomio entre el factor x t

Por ejemplo, el polinomio general:

En este punto, existen tres posibilidades:

El cociente es un polinomio de tercer grado o mayor. En tal caso, el mtodo de

1.3.1.1. PROGRAMA Y DIAGRAMA DE FLUJO