\documentclass[12pt]{article}
\usepackage{graphicx,amsmath,amsfonts,amssymb,epsfig,euscript}
\usepackage[T1]{fontenc}
\usepackage[utf8x]{inputenc}

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

\newcommand{\R}{\mathbb{R}}
\newcommand{\C}{\mathbb{C}}
\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}
\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\Huge {\sc Análisis Numérico}}}
 \vskip 0.2cm
 \centerline{\ttfamily Segundo Cuatrimestre 2016}
 \hrulefill

 \bigskip
 \centerline{\bf TP N$^\circ$ 1-B. Ecuación de Burgers, Método BDF con extrapolación}
 \bigskip

El objetivo de este trabajo es resolver la ecuación de Burgers en una dimensión
$$ u_t + u u_x = \mu u_{xx} $$
con condiciones de borde periódicas en $[0,2\pi]$, utilizando un método implícito en la difusión, y explícito en la parte no-lineal, de orden al menos 2 en el tiempo. 

Para la discretización espacial, se consideran los operadores usuales de diferencias centradas de orden 2, y 
para la discretización temporal buscamos desarrollar un método totalmente implícito (Ej. 11 Práctica 2), utilizando fórmulas de diferenciación 
a izquierda (o BDF) para aproximar $f'(t)$:
\begin{align*}
 \mbox{BDF1: } f'(t+h) \sim & \frac{1}{h} \left( f(t+h) - f(t) \right) \\
 \mbox{BDF2: } f'(t+h) \sim & \frac{1}{h} \left( \frac 32 f(t+h) - 2f(t) + \frac 12 f(t-h) \right) \\
 \mbox{BDF3: } f'(t+h) \sim & \frac{1}{h} \left( \frac{11}{6} f(t+h) - 3f(t) + \frac 32 f(t-h) - \frac 13 f(t-2h) \right) \\
 \mbox{BDF4: } f'(t+h) \sim & \frac{1}{h} \left( \frac{125}{60} f(t+h) - 4f(t) + 3f(t-h) - \frac{4}{3} f(t-2h) + \frac{1}{4} f(t-3h) \right)
\end{align*}
Así, se puede obtener un método totalmente implícito, de la forma 
$$\frac{\frac 32 u^{n+1}_i -2u^n_i + \frac 12 u^{n-1}_i}{\Delta t} = \delta_{xx}^2(u^{n+1}_i) + u^{n+1}_i\delta_{x}(u^{n+1}_i)$$
Para evaluar el término no-lineal a tiempo $n+1$ sin tener que resolver un sistema de ecuaciones no-lineales, se aproxima al mismo \textit{por extrapolación} utilizando $p$ valores anteriores -- 
la expresión que resulta de evaluar a tiempo $n+1$ el polinomio interpolador de Lagrange que utiliza los tiempos $n, n-1,\dots$ es:
$$ f(t+h) \sim \sum_{k=0}^{p-1} (-1)^k {{p}\choose{k+1}} f(t-kh), \qquad p\ge0 $$

\textbf{Para hacer (los primeros 3 items a modo de sugerencia):} 
\begin{enumerate}
 \item Implemente la fórmula de extrapolación con $p$ valores anteriores, y tome como \textit{test} un caso conocido. Comando útil: \texttt{nchoosek}.
 \item Implemente las fórmulas de derivación a orden 2, 3 y 4, y compare contra un caso conocido para comprobar su implementación.
 \item Implemente la rutina BDF2 para resover la ecuación del calor con condiciones de borde periódicas. 
 
 \textbf{Sugerencia:} Para implementar discretizaciones espaciales con borde periódico, 
 añada $1$s en las esquinas opuestas a la diagonal de la matriz usual que discretiza la derivada. 

 \item Combine las rutinas anteriores para obtener un esquema implícito en la difusión y explícito en la no-linealidad para la ecuación de Burgers. 
 Experimente con distintos valores de $s$ y $p$ para obtener un método estable y del orden más alto posible, y compare contra la solución exacta provista a continuación,
 para valores de $\mu << 1$.
\end{enumerate}

El problema  tiene como solución exacta
$$ u(x,t) = -\frac{2\mu}{\phi}{\phi_x} + 4, \qquad \mbox{con} \qquad \phi(x,t) = e^{\frac{-(x-4t)^2}{4\mu(t+1)}} + e^{\frac{-(x-4t-2\pi)^2}{4\mu(t+1)}} $$
si se toman como condiciones iniciales $u(x,0)$. 

\textbf{Sugerencia:} Para ahorrar tiempo en la implementación de esta función conocida, 
puede utilizar rutinas de cálculo simbólico como \texttt{syms,diff} y \texttt{matlabFunction}, como en el código de ejemplo provisto a continuación.
\begin{verbatim}
syms x nu t
phi = exp(-(x-4*t).^2/(4*nu*(t+1))) + exp(-(x-4*t-2*pi).^2/(4*nu*(t+1)));
phiprime = diff(phi,x);
usym = -2*nu/phi*phiprime+4;
sol_exacta = matlabFunction(usym); 
\end{verbatim}

 
 \textbf{(Opcional):} Un método alternativo para aplicar las discretizaciones espaciales con condiciones periódicas 
 es descomponer un vector $f$ en la base de Fourier discreta, aplicar el autovalor $-4\sin^2(\xi h)$ 
 a cada modo, y volver luego a la base original. Verifique que éste procedimiento coincide con aplicar la matriz $A$, y 
 luego obtenga un método de orden arbitrariamente alto, utilizando la derivada exacta 
 del modo de Fourier en vez del autovalor correspondiente al método de orden 2. 
 Los cambios de base pueden implementarse construyendo las matrices correspondientes, o bien utilizando las rutinas \texttt{fft,ifft} que requieren $O(n \log n)$ operaciones.
 


\end{document}