\documentclass[a4paper]{article}
\usepackage{amsmath,amsthm,amsfonts}
\usepackage{amssymb}
\usepackage{xcolor}
\usepackage[utf8]{inputenc}


\usepackage{listings}
\lstset{breaklines=true,tabsize=4,frame=single,numbers=left,,basicstyle=\ttfamily}
\definecolor{matlabgreen}{RGB}{34,139,34}
\lstset{language=matlab,
  numberstyle=\tiny\color{gray},
%  numbers=none,xleftmargin=10pt,
  stringstyle=\color{red},
  basicstyle=\ttfamily,
  keywordstyle=\color{blue},
  commentstyle=\color{matlabgreen},
  basewidth=0.45em,
  showstringspaces=false,
  captionpos=b,
  frame=single
}

\usepackage{graphicx}
\usepackage{subfigure}
\usepackage{caption}
\usepackage{euscript}
\usepackage{units}
\usepackage{enumitem}
%\usepackage[notref,notcite,color]{showkeys}
\usepackage[bookmarks=true,hidelinks]{hyperref}
\usepackage{cleveref}
\usepackage{bbm}
\usepackage{thmtools}
\usepackage{bm}
\usepackage{xcolor}
%\usepackage{calrsfs}
\usepackage{tikz}
\usepackage{mdframed}
\usepackage{titling}	% For titlingpage
\renewcommand*{\thefootnote}{\fnsymbol{footnote}}	% Use symbols for footnotes

\newtheorem*{theorem}{Theorem}
\theoremstyle{definition}
\newtheorem*{definition}{Definition}
\newtheorem*{remark}{Remark}
\newtheorem*{example}{Example}


\renewcommand{\geq}{\geqslant}
\renewcommand{\leq}{\leqslant}
\renewcommand{\epsilon}{\varepsilon}
\newcommand{\eps} {\varepsilon}
\renewcommand{\phi}{\varphi}
\newcommand{\R}{\mathbb{R}}
\newcommand{\Q}{\mathbb{Q}}
\newcommand{\N}{\mathbb{N}}
\newcommand{\Z}{\mathbb{Z}}
\renewcommand{\P}{\mathbb{P}}
\DeclareMathOperator{\sgn}{sgn}
\newcommand{\ind}{\mathbbm{1}} 
\newcommand{\ip}[2]{\langle#1,#2\rangle}
\newcommand{\trans}{{\mathsf{T}}}


\newcommand{\update}[1]{{\textcolor{red}{#1}}}


\newlist{problems}{enumerate}{3}
\setlist[problems,1]{label=\textbf{Problem \arabic*.},ref=\arabic*,align=left}
\setlist[problems,2]{label=(\alph*)}
\setlist[problems,3]{label=(\roman*)}

\newlist{problemsolutions}{enumerate}{1}
\setlist[problemsolutions,1]{label=(\alph*)}
	


\begin{document}


\begin{titlingpage}
    \vspace*{-8em}
    \setlength{\parskip}{1.5ex plus 0.5ex minus 0.2ex}
    \setlength{\parindent}{0pt}

    \begin{flushright}
        \small\today
    \end{flushright}

    \begin{center}
        \vspace{2em}
        {
            \bfseries\sffamily\huge
            MAT3110 -- Introduction to numerical analysis
            %% Alternatively:
            %% COURSEXXXX \textthreequartersemdash\ COURSE NAME
        }
        \vskip0.2ex
        \Large
        Mandatory assignment 2 of 2
        \vskip1ex
    \end{center}

    \subsubsection*{Submission deadline}
    Thursday 29 October 2020 at 14:30 in Canvas
    (\href{https://canvas.uio.no}{\underline{canvas.uio.no}}).

    \subsubsection*{Instructions}
    Bachelor students can choose between scanning handwritten notes or use \LaTeX. Master students are required to use \LaTeX.
    The assignment must be submitted as a single PDF file.
    Scanned pages must be clearly legible; please use either the ``Color'' or ``Grayscale'' settings, not ``Black and white''.
    The submission must contain your name, course and assignment number.

    It is expected that you give a clear presentation with all necessary explanations.
    Remember to include all relevant plots and figures.
Students who fail the assignment, but have shown that they have made a genuine attempt at solving the exercises, are given a second attempt at revising their answers.
    All aids, including collaboration, are allowed, but the submission must be written by you and reflect your understanding of the subject.
    If we doubt that you have understood the content you have handed in, we may request that you give an oral account.

    In exercises where you are asked to write a computer program, you need to hand in the code along with the rest of the assignment.

    \subsubsection*{Application for postponed delivery}
    If you need to apply for a postponement of the submission deadline due to illness or other reasons, you \textbf{have to} contact the student administration at the Department of Mathematics (e-mail: \href{mailto:studieinfo@math.uio.no}{studieinfo@math.uio.no}) well before the deadline.

    All mandatory assignments in this course must be approved in the same semester, before you are allowed to take the final examination.

    \subsubsection*{Complete guidelines for delivery of mandatory assignments:}
        \texttt{\href{/english/studies/examinations/compulsory-activities/mn-math-mandatory.html}
        {
            \resizebox{\textwidth}{!}
            {\underline{www.uio.no/english/studies/examinations/compulsory-activities/mn-math-mandatory.html}}
        }}
\end{titlingpage}

    

\noindent
\textbf{Grading:} To get your submission approved, you need to show that you have tried to solve every exercise, and have solved about $\unitfrac23$ of the problems correctly. Each of the three problems 1, 2, 3 are weighted by $\unitfrac13$ each.

\bigskip

In this assignment we will consider the \textit{ordinary differential equation}\footnote{This ODE can be viewed as a model of a piece of string. The ends of the string are fixed at the two points $(0,\alpha)$ and $(1,\beta)$, and $y(x)$ denotes the height above the ground of the string at the position $x$. The right-hand side $f$ is then a force acting in the vertical direction.}
\begin{equation}\label{eq:ode}
\begin{split}
y''(x) = f(y(x))& \qquad \text{for } x\in(0,1) \\
y(0)=\alpha,\ y(1)=\beta.&
\end{split}
\end{equation}
Here, $f\in C(\R)$ is a given continuous function, $\alpha,\beta\in\R$ are given ``boundary data'', and $y:[0,1]\to\R$ is the unknown function. A function $y=y(x)$ is a \textit{solution of \eqref{eq:ode}} if $y''(x)=f(y(x))$ at all points $x\in(0,1)$, and $y(0)=\alpha$, $y(1)=\beta$.

If $f$ is a linear or constant function then you have learnt how to solve \eqref{eq:ode} exactly. If $f$ is \textit{nonlinear} then a solution still exists, but in general it does not have a closed form expression. Our goal will be to numerically approximate this solution.

\begin{problems}[resume]
\item
Let $N\in\N$. We partition the domain $x\in[0,1]$ into $N$ intervals with endpoints
\[
x_i = ih \quad (i=0,1,\dots,N) \qquad \text{where } h = \tfrac{1}{N}.
\]
At each ``grid point'' $x_i$ we approximate $z_i \approx y(x_i)$.
\begin{problems}
\item
Explain why the set of equations
\begin{equation}\label{eq:dode}
\begin{split}
-z_{i-1}+2z_i-z_{i+1} = h^2f(z_i)& \qquad \text{for } i=1,2,\dots,N-1 \\
z_0=\alpha,\ z_N=\beta&
\end{split}
\end{equation}
gives an approximation of \eqref{eq:dode}.

\item
Denote the unknown by $z = \begin{pmatrix}z_0 & z_1 & \dots & z_N\end{pmatrix}^\trans$. Find a tridiagonal\footnote{A matrix $A$ is \textit{tridiagonal} if $A_{i,j}=0$ for $j>i+1$ and $j<i-1$.} matrix $A\in\R^{(N+1)\times(N+1)}$ and a vector $b(z)\in\R^{N+1}$ such that \eqref{eq:dode} can be equivalently written as
\begin{equation}\label{eq:dode2} %\tag{\ref{eq:dode}'}
Az = b(z).
\end{equation}
(Note that the vector $b(z)$ depends on $z$. Make sure that you include the boundary conditions from \eqref{eq:dode} in \eqref{eq:dode2}.)

\item 
A matrix $A\in\R^{n\times n}$ is \textit{diagonally dominant} if
\begin{enumerate}[label=(\roman*)]
\item none of the rows in $A$ contains only zeros
\item $|A_{i,i}| \geq \sum_{j\neq i} |A_{i,j}|$ for all $i=1,\dots,n$
\item there is at least one index $i$ such that $|A_{i,i}| > \sum_{j\neq i} |A_{i,j}|$. 
\end{enumerate}
It is a fact (which you do not need to prove) that all tridiagonal, diagonally dominant matrices are invertible. (The proof is an elementary, but tedious, proof showing that Gaussian elimination with $A$ works.)

Prove that $A$ from \eqref{eq:dode2} is invertible.

\item
Consider now a general linear system 
\begin{equation}\label{eq:linsys}
Bu=c, 
\end{equation}
where $c\in\R^n$ is given, $u\in\R^n$ is unknown, and where $B\in\R^{n\times n}$ is tridiagonal and diagonally dominant. Explain why solving \eqref{eq:linsys} using Gaussian elimination only requires $O(n)$ floating point operations (flops), as opposed to $O(n^3)$ for a general matrix $B$.
\end{problems}
\end{problems}



We now wish to solve \eqref{eq:dode2}. Since $b$ depends on $z$, the equation is nonlinear, and we cannot solve it directly. Instead, we will use an iterative solver which, given an initial guess $z^{(0)}\in\R^{N+1}$, computes approximations $z^{(1)},z^{(2)},\dots\in\R^{N+1}$. To measure the difference between successive iterations (which should be an indication of the error) we will use the $L^2$-type norm
\[
\big\|z^{(k+1)}-z^{(k)}\big\| = \sqrt{h\sum_{i=0}^N \big|z^{(k+1)}_i-z^{(k)}_i\big|^2}.
\]
We will use two iterative solvers: a fixed point iteration and Newton's method. 

\begin{problems}[resume]
\item\label{prob:fp}
Since $A$ is invertible, we can write \eqref{eq:dode2} as
\begin{equation}\label{eq:dode3}\tag{\ref{eq:dode2}'}
z = A^{-1}b(z).
\end{equation}

\begin{problems}
\item\label{prob:fpiter}
Use \eqref{eq:dode3} to find a fixed point iteration which approximates $z$.

\item 
Using a programming language of your choice, write a function\footnote{Functions can be passed as variables both in Matlab, Python, and other languages.} 
\begin{center}
\texttt{fixedpointODE(f, alpha, beta, N, tolerance)} 
\end{center}
which takes as input $f$, $\alpha$, $\beta$, $N$, and a positive number \texttt{tolerance}, runs the fixed point iteration from problem \ref{prob:fpiter} until $\big\|z^{(k)}-z^{(k-1)}\big\|<\mathtt{tolerance}$, and returns the final iterate $z^{(k)}$.

\item\label{prob:fptest}
Test your code on the two problems
\begin{subequations}
\begin{align}
&\alpha=0,\ \beta=1,\ f(z)=-1\label{eq:prob1} \\
\intertext{and}
&\alpha=1,\ \beta=1,\ f(y)=-y^2,\label{eq:prob2}
\end{align}
\end{subequations}
in both cases using $N=50$, $z^{(0)}=0$ and $\mathtt{tolerance}=10^{-12}$. Include plots of $z$ as a function of $x$ from both runs.  (For \eqref{eq:prob1} you can compute the exact solution by hand, in case you want to verify your code.)

\item Estimate how many floating point operations (flops) each iteration requires. You only need a rough estimate, such as $O(N^p)$ for some $p$.
\end{problems}



\item\label{prob:newton}
Alternatively, we can write \eqref{eq:dode2} as
\begin{equation}\label{eq:dode4}\tag{\ref{eq:dode2}{'}'}
g(z)=0, \qquad \text{where } g(z)=Az-b(z).
\end{equation}
\begin{problems}
\item \label{prob:newtoniter}
Write down Newton's method for solving \eqref{eq:dode4}.

\item 
Write a function 
\begin{center}
\texttt{newtonODE(f, df, alpha, beta, N, tolerance)} 
\end{center}
which takes as input $f$, its derivative $f'$, as well as $\alpha$, $\beta$, $N$, and a positive number \texttt{tolerance}, runs the Newton iteration from problem \ref{prob:newtoniter} until $\big\|z^{(k)}-z^{(k-1)}\big\|<\mathtt{tolerance}$, and returns the final iterate $z^{(k)}$.

\item Test your code in the same way as in problem \ref{prob:fp}\ref{prob:fptest}.

\item Estimate how many floating point operations (flops) each iteration requires. You only need a rough estimate, such as $O(N^p)$ for some $p$.

\item Based on how many iterations your test runs required, which of the fixed point iteration or Newton's method do you think is faster?
\end{problems}
\end{problems}



\end{document}
