\documentclass{article}
\usepackage{graphicx}  % To manage external pictures
\usepackage{float}
\usepackage[colorlinks=true,linkcolor=blue,citecolor=blue,urlcolor=blue]{hyperref}
\usepackage[left=3cm,top=3cm,bottom=3.5cm,right=3cm]{geometry}	% For easy management of document margins
\parindent 0em

\begin{document}

<<echo=FALSE>>=
cat("UNIVERSIDAD NACIONAL AGRARIA\nDPTO.Estadistica e Informatica\nMetodos Numericos y Simulacion")
@
$$Final,~ Diciembre ~del ~2014$$

\begin{enumerate}

\item Dada la matriz A sim\'etrica, Aplique el m\'etodo de Jacobi para hallar las ra\'ices y vectores caracter\'isticos, si:

$$A=\left (\begin{array}{ccc} 4&0&0 \\ 0&5&2 \\ 0&2&1 \end{array} \right )$$

por Jacobi, la matriz de rotaci\'on plana es:

$$P=\left (\begin{array}{ccc} 1&0&0 \\ 0&C&-S \\ 0&S&C \end{array} \right )$$

C=coseno
S=seno

El \'angulo $\theta$ de la rotaci\'on plana es:

$\theta=0.5~atan\left(\frac{2a_{23}}{a_{22}-a_{33}}\right)$

<<echo=false>>=
A<-rbind(c(4,0,0),c(0,5,2),c(0,2,1))
# Hacer cero a(2,3)=a(3,2)=0, i=2,j=3
ang<-0.5*atan(2*A[2,3]/(A[2,2]-A[3,3]))
C<-cos(ang)
S<-sin(ang)
P<-rbind(c(1,0,0),c(0,C,-S),c(0,S,C))
A1<-t(P)%*%A%*%P
A1
# valores caracteristicos
round(diag(A1),4)
# vectores caracteristicos
P
@

Esto significa que para $\lambda_1 = 4$ el vector asociado es $\phi_1=\left( \begin{array}{c}1\\0\\0 \end{array}\right)$

Para $\lambda_2 = 5.8284$, el vector asociado $\phi_2=\left( \begin{array}{c}0\\0.9238\\0.3826 \end{array}\right)$

y para $\lambda_3 =0.1715$, el vector asociado $\phi_3=\left( \begin{array}{c}0\\-0.3826\\0.9238 \end{array}\right)$ 

\----

\item Encontrar la relaci\'on de las diferencias dividas con diferencias finitas hacia adelante si las $X's$ est\'an igualmente espaciadas. Utilice n pares de puntos $(x_i, f_i)$

Para diferencias finitas hacia adelante, los puntos $x_i$ est\'an igualmente espaciadas  $x_{k+1} = x_k + h$, a las que corresponde valores de funci\'on $f(x): f(x_k)$ y se definen las primeras diferencias finitas hacia delante como: 

$\Delta f_k=f_{k+1}-f_k$

Para el caso de las diferencias divididas, los espaciamientos pueden ser diferentes, cuando los espaciamientos son iguales se puede establecer una equivalencia con las diferencias finitas hacia adelante, as\'i:

$\{x_0,x_1\}=\frac{f(x_0)-f(x_1)}{x_0-x_1}=\frac{-\Delta f_0}{-h}=\frac{\Delta f_0}{h}$

$\{x_0,x_1,x_2\}= \frac{\{x_0,x_1\}-\{x_1,x_2\}}{x_0-x_2}=\frac{\frac{\Delta f_0}{h}-\frac{\Delta f_1}{h}}{2h}=\frac{\Delta^2 f_0}{2h^2}$

$\{x_0,x_1,x_2,x3\}= \frac{\{x_0,x_1,x2\}-\{x_1,x_2,x_3\}}{x_0-x_3}=\frac{\frac{\Delta~^2 f_0}{2h^2}-\frac{\Delta~^2 f_1}{2h^2}}{3h}=\frac{\Delta~^3 f_0}{3*2h^3}=\frac{\Delta~^3 f_0}{3!h^3}$

Se supone cierto para $(n-1)$

$\{x_0,x_1,x_2,...x_{n-1}\}=\frac{\Delta~^{n-1} f_0}{(n-1)!h^{n-1}}$

Entonces para $n$ puntos

$\{x_0,x_1,x_2,...,x_n\}= \frac{\{x_0,x_1,x_2,...,x_{n-1}\}-\{x_1,x_2,x_3,...,x_n\}}{x_0-x_n}$

Reemplazando sus equivalentes:

$\{x_0,x_1,x_2,...x_n\}= \frac{\frac{\Delta~^{n-1} f_0}{(n-1)!h^{n-1}} -\frac{\Delta~^{n-1} f_1}{(n-1)!h^{n-1}}}{nh}=\frac{\Delta~^n f_0}{n!h^n}$

Generalizando se tiene:

$\{x_i,x_{i+1},...x_{i+n}\}=\frac{\triangle^n f_i}{n! h^n}$

\----

\item Hallar las diferencias divididas hasta segundo orden de la funci\'on $e^{-x^2}$ si se tabula la funci\'on para $x=c(0.1,~ 0.2,~ 1,~ 1.5)$. 

Realice la interpolaci\'on para $x=0.5$.

$f(x)=e^{-x^2}$

La funci\'on tabulada para $x=c(0.1,~ 0.2,~ 1,~ 1.5)$ es:

<<echo=FALSE>>=
f<-function(x)exp(-x^2)
x<-c(0.1, 0.2, 1, 1.5)
cbind(x,f(x))
@

el verdadero valor es: $f(0.5)=0.7788008$

Aplicando diferencias divididas:

<<echo=FALSE>>=
y<-NULL # primera diferencia
z<-NULL # segunda diferencia
for(i in 1:3) y[i]<-(f(x[i+1])-f(x[i]))/(x[i+1]-x[i])
for(i in 1:2) z[i]<-(y[i+1]-y[i])/(x[i+2]-x[i])
# interpolaci\'on para x=0.5
w1<-f(x[1])+y[1]*(0.5-x[1])+z[1]*(0.5-x[1])*(0.5-x[2])
w2<-f(x[2])+y[2]*(0.5-x[2])+z[2]*(0.5-x[2])*(0.5-x[3])
cat("Primeras diferencias dividias\n")
print(y)
cat("Segundas diferencias dividias")
print(z)
cat("Valor interpolado para x=0.5")
cat("Iniciando en x0=",x[1],"es:",w1)
cat("Iniciando en x0=",x[2],"es:",w2)
@

Entonces el error absoluto es:

$\epsilon_{f(0.5)}=|0.8132038-0.7788008| = 0.0344$

$\epsilon_{f(0.5)}=|0.7135047-0.7788008| = 0.0652$

\--- 

\item Se desea hallar la integral de la funci\'on coseno entre $\pi/4$ y $\pi/2$. Utilice Simpson con extrapolaci\'on. Utilice el espaciamiento que mejor se adecue para su prop\'osito y determine las cifras significativas de esta aproximaci\'on respecto al verdadero valor.

\begin{figure}
\begin{center}
<<f1,echo=FALSE,fig=TRUE,height=3,width=6>>=
par(mar=c(2,2,0,0),cex=1)
f<-function(x) cos(x)
intervalo<-c(pi/4,pi/2)
curve(f,intervalo,col="blue",xlim=c(0,pi/2),type="h",axes=FALSE)
points(intervalo,f(intervalo),type="h",col="blue",lty=4,lwd=2)
I<-integrate(f,pi/4,pi/2)$value
Gintegra<-expression(integral(~~cos(x)~~dx,pi/4,pi/2) )
text(0.4,0.5,Gintegra)
text(0.4,0.3,substitute(paste("=",Iv), list(Iv=I)))
x0<-expression(pi/4)
xn<-expression(pi/2)
axis(1,c(pi/4,pi/2),c(x0,xn))
axis(2)
@
\caption{Integral\label{fig:f1}}
\end{center}
\end{figure}

En Simpson se requiere un n\'umero m\'inimo de 3 puntos para el c\'alculo del \'area, entonces se debe generar 5 puntos $(x_1,x_2,x_3,x_4,x_5)$, igualmente espaciados, de tal forma que se pueda utilizar $x_1,x_3,x_5$ para el \'area con $2h$ y los puntos  $(x_1,x_2,x_3)$ y $(x_3,x_4,x_5)$ para el \'area con espaciamientos iguales a $h$.

el valor de h seria:

$h=\frac{\frac{1}{2}\pi-\frac{1}{4}\pi}{4}=\frac{1}{16}\pi$

Los puntos de x, desde $pi/4$ hasta $pi/2$ serian:

$x_1=\frac{4}{16}\pi;~~$
$x_2=\frac{5}{16}\pi;~~$
$x_3=\frac{6}{16}\pi;~~$
$x_4=\frac{7}{16}\pi;~~$
$x_5=\frac{8}{16}\pi$

Espaciamiento=h
$I1=\frac {h}{3} \left(cos(x_1)+4cos(x_2)+2cos(x_3)+4cos(x_4)+cos(x_5) \right)$

Espaciamiento=2h
$I2=\frac {2h}{3} \left(cos(x_1)+4cos(x_3)+cos(x_5)\right)$

Extrapolaci\'on de Richardson
$area=\frac {2^4~I1-I2}{2^4-1}$

N\'umero de cifras significativas

Si $I$ es el valor exacto de la integral, y $area$ el aproximado, entonces:

$\epsilon_area=|area-I|$

Aplicando la relaci\'on de cifras significativas se tiene:

$\epsilon_ area \le 0.5*10^{m-n+1}$

Como $m=-1$ para este caso, entonces:

$\epsilon_area \le 0.5*10^{-n}$

<<echo=FALSE>>=
par(mar=c(2,2,1,0),mfrow=c(1,2),cex=1)
x<-seq(pi/4,pi/2,length=5)
y<-cos(x)
#cbind(x,y)
f<-function(x) cos(x)
h<-(pi/2-pi/4)/4
# para 2h 
I2<-2*h/3*(f(x[1])+4*f(x[3])+f(x[5]))
# para h
I1<-h/3*(f(x[1])+4*f(x[2])+2*f(x[3])+4*f(x[4])+f(x[5]))
cat("I1:",I1,"\n")
cat("I2:",I2,"\n")
area<-(2^4*I1-I2)/(2^4-1)
cat("Area:",area,"\n")
cat("Error absoluto:",abs(area-I),"\n")
@

Como el error es menor de: $0.37~10^{-7}$; $0.37 < 0.5$;  implica que $10^{-7} = 10^{-n}$

Entonces $n=7$; por lo tanto la aproximaci\'on tiene 7 cifras significativas.

\---

\item Resolver el siguiente sistema de ecuaciones diferenciales por el m\'etodo de Euler:

$$Y''' - 2XY''+ YY'= e^x$$

si las condiciones iniciales son:

$Y(1)=0.5$

$Y'(1)=-0.5$

$Y''(1)=0.25$

Hallar $Y(1.1),Y(1.2),Y(1.3)$, utilice $h=0.1$

Este sistema se debe transformar en un sistema de ecuaciones de primer orden.

el sistema original es:

$$Y'''= e^x + 2XY''- YY'= f(x,Y,Y',Y'')$$

la transformaci\'on ser\'a:

$Z1 = Y$

$Z2 = Y'$

$Z3 = Y''$

Entonces:

$Z1'=Z2$

$Z2'=Z3$

$Z3'= f(x,Z1,Z2,Z3)$

valores iniciales:

$Z1(1)=0.5$; $Z2(1)=-0.5$; $Z3(1)=0.25$; 

el proceso iterativo segun Euler para h=0.1

$Z1_{i+1}=Z1_i+0.1~Z2_i$

$Z2_{i+1}=Z2_i+0.1~Z3_i$

$Z3_{i+1}=Z3_i+0.1~(e^x_i+2x~Z3- Z1~Z2)$

Mediante R.

<<echo=FALSE>>=
x<-NULL;Z1<-NULL; Z2<-NULL; Z3<-NULL
# valores inciales
x[1]<-1; h=0.1
Z1[1]<-0.5; Z2[1]=-0.5;Z3[1]=0.25
for(i in 1:3) {
Z1[i+1]<-Z1[i]+0.1*Z2[i]
Z2[i+1]<-Z2[i]+0.1*Z3[i]
Z3[i+1]<-Z3[i]+0.1*(exp(x[i])+2*x[i]*Z3[i]-Z1[i]*Z2[i])
x[i+1]<-x[i]+h
}
cbind(x,Z1,Z2,Z3)
# La solucion corresponde a Y=Z1
cbind(x,Y=Z1)
@

\end{enumerate}

Puntaje 4 c/u
\end{document}