\documentclass[12pt]{article}
%\include{amsfonts}
\usepackage{amssymb,amsmath,latexsym,epsfig,euscript,multicol}
\usepackage{enumitem}
% \usepackage[utf8x]{inputenc}
\usepackage[latin1]{inputenc}

% Caracteres especiales
\def\A{\mathbb{A}}
\def\C{\mathbb{C}}
\def \N{\mathbb{N}}
\def \P{\mathbb{P}}
\def \Q{\mathbb{Q}}
\def \R{\mathbb{R}}
\def \Z{\mathbb{Z}}
\def \sen{\textrm{sen}}

\newtheorem{ejer}{Ejercicio}
\newcommand{\bej}{\begin{ejer}\rm}
\newcommand{\fej}{\end{ejer}}

\def\dt{\Delta t}
\def\dx{\Delta x}

\topmargin-2cm \vsize 29.5cm \hsize 21cm
\setlength{\textwidth}{16.75cm}\setlength{\textheight}{23.5cm}
\setlength{\oddsidemargin}{0.0cm}
\setlength{\evensidemargin}{0.0cm}


\begin{document}

%----------------------------------------------------------
%Encabezado:

\centerline{{\small Universidad de Buenos Aires - Facultad de Ciencias Exactas y Naturales - Depto. de Matemática}}
 
 \vskip 0.2cm
 \hrulefill
 \vskip 0.2cm

 \centerline{{\bf\Large{\sc Elementos de Cálculo Numérico / Cálculo Numérico}}}
 \vskip 0.2cm
 \centerline{\ttfamily Curso de Verano 2019}
 \hrulefill

 \bigskip
 \centerline{\bf Práctica N$^\circ$ 2: Ecuaciones Diferenciales: Problemas de valores iniciales.}
 \bigskip


%----------------------------------------------------------
% 1
\bej
 Se considera la siguiente ecuación diferencial:
$$
\left\{
\begin{array}{rl}
y'(t) = & 2y(t) - 5\sen(t)\\
y(0) =  & 1
\end{array}
\right.
$$
cuya solución exacta es la función $y(t) = 2\sen(t) +
\cos(t)$. 
\begin{enumerate}[label=\alph*)]
\item Escribir la iteración del método de Euler para esta ecuación. 
\item Calcular el error de truncado local. 
\item ¿Qué paso $h$ debe elegirse para que el error al estimar $y\big(\frac{\pi}{2}\big)$ sea menor que $10^{-2}$?
\end{enumerate}
\fej

%----------------------------------------------------------
% 2
\bej Graficar simultáneamente en la región $[0, 10]\times [0, 10]$ las soluciones que se obtienen del 
problema de valores iniciales
$$\begin{array}{rcl} y'(t) & = & (y(t) - 5).(\cos^2(t) - 0.5),\\y(0) & = & k,\end{array}$$
al utilizar el método de Euler con paso $h = 0.01$ para $k = 0, 1, \dots, 10$.
\fej

%----------------------------------------------------------
% 3
\bej Considerar el problema $y'=-2ty$, $y(0)=1$, con $t\ge 0$.
\begin{enumerate}[label=\alph*)]
\item Determinar una cota, en términos de $h$, para el error
cometido si se usa el método de Euler para calcular $y(1)$.

\item ¿Cómo deberí­a tomar $h$ si se desea que el error
cometido sea menor que $10^{-2}$?

\item Calcular la solución en $t=1$ usando el valor de $h$
obtenido en el item previo, y verificar las estimaciones previstas comparando con la solución exacta.
\end{enumerate}
\fej

%----------------------------------------------------------
% 4
\bej \label{ej:eulersimple}  Escribir un programa que implemente el método de Euler explí­cito para resolver ecuaciones de la forma 
$$\begin{array}{rcl} y' & = & f(t,y)\\ y(t_0) & = & y_0,\end{array}$$ tomando como parámetros la función $f$, los tiempos inicial y final $t_0$ y $t_f$, el paso $h$ y el dato inicial $y_0$; y arrojando como resultados el vector $t=(t_0,t_0+h,\dots,t_f)$ y la solución $y$. 
\fej

%----------------------------------------------------------
% 5
\bej Se quiere verificar numéricamente el orden de convergencia de los métodos de Euler y Taylor de orden 2. 
Para ello: resolver numéricamente el problema $y'=y, \, y(0) =1$, en el intervalo $[0,1]$ con ambos métodos, 
tomando $h=0.1$, $h=0.0625$, $h=0.05$, $h=0.025$ y $h=0.01$. Para cada $h$ calcular el error que se comete al 
aproximar $y(1)$: $e_h = |y(1)-y_N|$. Graficar $\log(e_h)$ en función de $\log(h)$. ¿Qué se espera ver? ¿El resultado es consistente con el esperado?
\fej

%----------------------------------------------------------
% 6
\bej Considerar el problema $\left\{\begin{array}{ll}
y'=\lambda y \\
y(0)=y_0 \end{array}\right.$.
 
\begin{enumerate}[label=\alph*)]
\item Verificar que el método de Euler con paso $h$ genera la
sucesión:
\[ y_i=(1+\lambda h)^iy_0 \ \ i=0,1,\ldots.  \]
\item  Para $\lambda <0$, determinar para qué valores de $h$
ocurre que $y_i \rightarrow 0$ cuando $i \to \infty$. Comparar con la solución exacta. 
\item Resolver usando el programa del Ejercicio \ref{ej:eulersimple} para distintos valores de $\lambda$ ($\lambda=1,10,50,100$) y comparar con la solución exacta. ¿Qué sucede?
\item Repetir los items anteriores considerando el método de Euler implí­cito. ¿Qué se observa?
\end{enumerate}
\fej

%----------------------------------------------------------
% 7
\bej Considerar la ecuación: 
$$\left\{\begin{array}{rcl} y'(t) & = &\dfrac{e^{-y}}{t}\\ y(1) &= &0.\end{array}\right.$$
\begin{enumerate}[label=\alph*)]
\item Probar que $0\le y(t)\le t$ para $t\ge 1$. 
\item Escribir la iteración dada para esta ecuación por el método de Euler. Probar que la solución numérica resultará creciente. 
\item Calcular el error de truncado del método de Euler aplicado a la ecuación.
\item Dar un valor de paso $h$ que garantice que el error de la estimación numérica de $y(2)$ sea menor que $10^{-3}$.
\end{enumerate}
\fej

%----------------------------------------------------------
% 8
\bej Probar que una ecuación de orden $n$:
$$y^{(n)} = f\big(t,y,y',\dots,y^{(n-1)}\big),$$ 
se puede escribir como un sistema de $n$ ecuaciones de primer orden. Mostrar que un problema de valores iniciales para la primera se transforma en un problema de valores iniciales para el sistema.
\fej

%----------------------------------------------------------
% 9
\bej
Modificar el programa del Ejercicio \ref{ej:eulersimple} para que acepte ecuaciones vectoriales: la solución $y$ deberá ser una matriz de $m\times n$, donde $m$ es el número de pasos temporales y $n$ la cantidad de variables del problema. De este modo, la fila $i$ de $y$ corresponderá al valor de la solución en todas sus variables a tiempo $t_i$ 
\fej

%----------------------------------------------------------
% 10
\bej{\bf Sistema predador-presa:} \label{ej:predador-presa}
Se tienen dos poblaciones, una de predadores y otra de presas, cuyo número a tiempo $t$ denotamos $x(t)$ e $y(t)$ respectivamente. En ausencia de presas, $x$ tiende a decaer a una tasa $\alpha$, mientras que en ausencia de predadores $y$ tiende a crecer a una tasa $\beta$. Además, los encuentros de predadores y presas hacen crecer la población de los primeros y decrecer la de los segundos, de acuerdo a cierta proporción. De esto modo, se obtiene el sistema:
$$\begin{array}{rcl}
\dot{x} &=& -\alpha x + \gamma xy \\
\dot{y} &=& \beta y - \delta xy,
\end{array}$$
donde $\gamma y$ es la tasa de crecimiento de $x$ (mayor cuanto más presas haya) y $\delta x$ es la tasa de mortandad de presas (mayor cuanto más predadores haya). Se asume que $\alpha,\beta,\gamma,\delta$ son todos positivos. 

\begin{enumerate}[label=\alph*)]
\item Dar condiciones sobre los parámetros y los niveles de $x$ e $y$ que garanticen la estabilidad de las poblaciones. Es decir, que $x(t+\Delta t)=x(t)$ e $y(t+\Delta t) =y(t)$ para todo $\Delta t>0$. 
\item Elegir valores de $\alpha$, $\beta$, $\gamma$, $\delta$, $x_0$ e $y_0$ que satisfagan las condiciones del item anterior y resolver utilizando el método de Euler. Realizar dos gráficos: uno de $x$ e $y$ en función de $t$ (simultáneamente) y otro de $y$ en función de $x$. 
¿Se mantiene constante la solución? 
\item Tomar $\alpha = 0.25$, $\beta=1$, $\gamma=\delta=0.01$, $x_0 = 80$ e $y_0=30$. Resolver utilizando el método de Euler y realizar gráficos como los del ítem anterior. 
\end{enumerate}
\fej

%----------------------------------------------------------
% 11
\bej
Considerar el método de Euler modificado:
$$y_{i+1} = y_i + h~f\left(t_i+\frac{h}{2},y_i+\frac{h}{2}f(t_i,y_i)\right).$$
Probar que el error de truncado es $O(h^2)$. ¿Qué ventaja presenta este método respecto del método de Taylor de segundo orden?
\fej

%----------------------------------------------------------
% 12
\bej{\bf Galileo:} Leer el siguiente párrafo:
\begin{quotation}
-- Pero, Simplicio, tengo la esperanza de que no seguirás el ejemplo de muchos otros que desvían la discusión de un punto principal y dicen que algunas de mis afirmaciones se apartan de la verdad por un cabello, y por este cabello esconden las faltas de otras teorías tan gruesas como un cable de navío. Aristóteles dice que `una esfera de hierro de 100 libras, cayendo desde una altura de 100 codos, llega a la tierra antes que una bola de 1 libra haya caído un simple codo'. Yo digo que las dos llegan al mismo tiempo. Tú encuentras, al hacer la experiencia, que la más pesada adelanta a la más ligera en 2 ó 3 dedos; ahora, no puedes esconder detrás de estos dos dedos los 99 codos de Aristóteles, ni puedes mencionar mi error y, al mismo tiempo, pasar en silencio el tuyo, mucho mayor.
\begin{flushright}
 Salviati, en \emph{Diálogo sobre dos nuevas ciencias} - Galileo Galilei.
 \end{flushright}
\end{quotation}

Viviani, estudiante de Galileo, afirma que su maestro realizó efectivamente el experimento descripto en el párrafo anterior, arrojando desde lo alto de la torre de Pisa una bala de cañón y una bala de mosquete. El objetivo de este ejercicio es reproducir numéricamente la experiencia de Galileo.

La posición de un objeto en caída libre puede modelarse con la ecuación:
\begin{equation}\label{eq}
 \ddot{x} = \frac{\gamma}{m}\dot{x}^2-g
\end{equation}
siendo $x$ la altura, $m$ la masa del cuerpo, $g=9.81 \frac{\rm{m}}{\rm{s}^2}$ la aceleración gravitatoria y $\gamma$ una constante que representa el rozamiento con el fluido en que se produce la caída. 
Deben darse condiciones sobre la altura y la velocidad iniciales. 

La Torre de Pisa mide $55.8$ mts. La masa de una bala de cañón es de $16$ Kg, y la de una bala de mosquete $0.0082$ Kg. Las costantes de rozamiento para cada bala son: $\gamma_c = 0.0058$ y $\gamma_m = 3.74\times 10^-5$, respectivamente (la diferencia se debe a la diferencia de tamaños).

Implementar un programa llamado \verb+galileo+ para obtener la dinámica de la caída de ambas balas utilizando el método de Euler modificado, y 
graficar, en una misma figura, la posición de cada bala en función del tiempo. A partir de los resultados obtenidos, responder:
\begin{enumerate}[label=\alph*)]
 \item ¿Cuánto tiempo tarda cada bala en tocar el suelo?
 \item Modificar el programa para que se detenga en el momento en que la bala cañón alcanza el suelo. ¿Cuán lejos del piso está la bala de mosquete?

 \emph{Nota:} No debe cometerse el mismo error que Simplicio al juzgar los resultados. La bala de cañón es alrededor de 2000 veces más pesada que la de mosquete. Consecuentemente, Aristóteles hubiese pronosticado que al llegar la bala de cañón al piso, la de mosquete habría descendido apenas 2 cm.
\end{enumerate}
\fej

%----------------------------------------------------------
% 13
\bej\label{ej:RK4}
Implementar un programa que resuelva sistemas de la forma:
$$y' = f(t,y),$$ utilizando el método de Runge Kutta de orden 4 dado por:
$$y_{i+1} = y_i + \frac{h}{6}\Big[K_1 +2K_2 +2K_3 +K_4\Big],$$
donde:
$$
\begin{array}{c}
K_1 = f(t_i,y_i),\\
K_2 = f\big(t_i+\frac{h}{2},y_i+\frac{h}{2}K_1\big),\\
K_3 = f\big(t_i+\frac{h}{2},y_i+\frac{h}{2}K_2\big),\\
K_4 = f\big(t_i+h,y_i+hK_3\big).
\end{array}
$$
Utilizar este método para resolver nuevamente el Ejercicio \ref{ej:predador-presa}. Comparar la solución con la obtenida con el método de Euler.
\fej


%----------------------------------------------------------
% 14
\bej\label{ej:tirooblicuo}{\bf Tiro oblicuo:}
Un proyectil de masa $m$ se arroja desde un punto del plano $(x_0,y_0)$, con una velocidad inicial dada por el vector $(v^x_0,v^y_0)$. La trayectoria del proyectil se rige por las ecuaciones dadas por la segunda ley de Newton:
$$\begin{array}{rcl}
m\ddot{x} &=& -\gamma \dot{x}\\ 
m\ddot{y} &=& -mg -\gamma \dot{y},
\end{array}$$
donde $g$ es la aceleración gravitatoria $g=9.81\frac{\textrm{m}}{\textrm{s}^2}$, y $\gamma$ es una constante de rozamiento con el medio en que se realiza el lanzamiento. 
Formular el problema en forma de sistema de orden uno. 

Tomando $m=10 \textrm{Kg}$ y $\gamma=0.2 \frac{\textrm{Kg}}{\textrm{s}}$, y suponiendo que el proyectil se lanza desde $30$ metros de altura con una velocidad inicial horizontal de $40\frac{\textrm{m}}{\textrm{s}}$, ¿qué distancia recorre antes de tocar el piso?

Hacer un programa que permita responder esta pregunta, utilizando el método de Euler modificado para resolver el sistema. 
\fej

%----------------------------------------------------------
% 15
\bej\label{pendulo}{\bf Péndulo:} Se considera el problema del péndulo:
$$\left\{\begin{array}{ll}
\ddot{\theta}(t)=-A\sen( \theta(t)) \mbox{ en } [0,T]\\
\dot{\theta}(0) = v_0 \\
\theta(0)=\theta_0 
\end{array}\right.$$
donde $\theta$ representa el ángulo que forma la vara del péndulo con la vertical.

\begin{enumerate}[label=\alph*)]
\item Formular el problema como un sistema de ecuaciones de orden uno.
\item Utilizar el método de Euler modificado, con paso $h=0.05$ para obtener una aproximación de la soluci\'on y graficarla.  
\item Graficar la solución que se obtiene al utilizar método de Runge Kutta del Ejercicio \ref{ej:RK4}.
\end{enumerate}

Pueden considerarse, a modo de ejemplo, los valores $A=7$, $T=10$, $\theta_0=\frac{\pi}{4}$, $v_o=0$.
\fej


\end{document}
