\documentclass[11pt,norsk,a4paper]{article}
\usepackage[T1]{fontenc}
\usepackage[utf8]{inputenc}
\usepackage{babel}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{amsthm}
\usepackage{subfigure}
\usepackage{../MAT-INF1100-oppg}
\usepackage{xcolor}
\usepackage{tikz}
\usetikzlibrary{decorations.pathreplacing}
\usepackage[hidelinks]{hyperref}
%\usepackage{wrisym}
%\usepackage{lucidabr}
\textwidth=14cm \oddsidemargin=1cm \evensidemargin=0cm		% A4
%\newif\ifpdf
%\ifx\pdfoutput\undefined
%\pdffalse % we are not running PDFLaTeX
%\else
%\pdfoutput=1 % we are running PDFLaTeX
%\pdftrue
%\fi
%
%\ifpdf
%\usepackage[pdftex]{graphicx}
%\else
%\usepackage{graphicx}
%\fi

\usepackage{listings}

\lstdefinelanguage{matlab}
{
  morekeywords={if,else,end,for,while,do,def,break,return},
  sensitive=false,
  morecomment=[l]{//}
}
  
\lstset{language=matlab,
  backgroundcolor=\color[rgb]{.95,.95,.95},
  numbers=none,xleftmargin=10pt,
  numberstyle=\tiny,stepnumber=5,numbersep=5pt,
  stringstyle=\color{red},
  basicstyle=\ttfamily,
  keywordstyle=\color{blue},
  commentstyle=\color{green},
  basewidth=0.45em,
  showstringspaces=false,
  captionpos=b,
  frame=single
}

\title{Obligatorisk oppgave2 i MAT-INF 1100 h酶sten 2017}


\begin{document}
%\begin{flushright} \small 1.~november, 2017
%     \end{flushright}
%\vspace{10mm}

\begin{center}\huge MAT-INF 1100: Obligatorisk oppgave 2\\
     \vspace{5mm}
\Large  Innleveringsfrist: torsdag 8. november 2018 kl. 14:30
\end{center}

\section*{}
Obligatoriske oppgaver (芦obliger禄) er en sentral del av MAT-INF1100 og er utmerket trening i 氓 besvare en matematisk oppgave.
Denne veiledningen er ment for 氓 klargj酶re noen krav vi som retter obligene stiller, b氓de for 氓 forenkle v氓r jobb i retteprosessen og for at dere skal slippe 氓 m氓tte levere p氓 nytt p氓 grunn av formaliteter.

For 氓 f氓 godkjent en oblig b酶r man ha f氓tt til om lag 2/3 av oppgavesettet og gitt alle deloppgavene et seri酶st fors酶k. 
Dersom det allikevel skulle v忙re en oppgave du overhodet ikke f氓r til, pr酶v 氓 sette ord p氓 hva du har tenkt og hvor det ikke lot seg l酶se.
Dersom man f氓r 芦ikke godkjent禄 p氓 en oblig og f氓r et nytt fors酶k vil det ha innleveringsfrist to uker etter f酶rste innlevering.  

\section{Krav til innlevering}

\begin{itemize}
	\item Besvarelsen skal leveres elektronisk via Devilry (\url{https://devilry.ifi.uio.no/}). Du velger selv om du skriver besvarelsen for h氓nd (og scanner besvarelsen) eller om du skriver l酶sningen direkte inn p氓 datamaskin (f.eks. ved bruk av \LaTeX). Skannede ark m氓 v忙re godt lesbare. Godkjente filtyper: PDF

	\item Hoveddelen av besvarelsen skal leveres som 茅n fil og i PDF-format. Det er lov 氓 skanne inn h氓ndskrevne ark s氓 lenge disse er godt lesbare. 

	\item Besvarelsen skal inneholde navn, emne og oblignummer.

	\item Husk at andre skal lese og forst氓 din besvarelse! Besvarelsene b酶r derfor f酶res p氓 en oversiktlig m氓te i den rekkef酶lgen de st氓r oppgitt.

	\item Alle plott og figurer som er med p氓 氓 besvare oppgaven m氓 inkluderes (husk ogs氓 riktige akser og enheter). Du b酶r ogs氓 skrive kort hva hver figur viser, for eksempel:
	\textit{Figur 2: Plott av funksjonen $f(x) = x^2 -3$.}

	\item Kode kan leveres ved siden av besvarelsen, men \textbf{all diskusjon} og \textbf{all besvarelse} av oppgaven skal inn i hovedfila. 
	Det vil f.eks. si at du ikke kan levere et Python-skript som besvarelse p氓 en oppgave hvor du blir bedt om 氓 lage en algoritme. Det er viktig at programkoden du leverer inneholder et kj酶reksempel, slik at det er lett 氓 se hvilket resultat programmet gir. 

	\item Det oppmuntres til 欧洲杯在线买球_欧洲杯投注网站推荐@ om oppgavene, men du skal selv ha formulert og skrevet den besvarelsen som leveres inn, og den skal gjenspeile din forst氓else av stoffet. Du kan bli bedt om 氓 redegj酶re muntlig for innholdet i din besvarelse.

	\item Besvarelsen b酶r inneholde diskusjon om resultatene dere kommer fram til. Eksempler p氓 sp酶rsm氓l man kan stille seg for 氓 vite hva man skal kommentere er: \textit{Gir dette plottet noen mening?} \textit{Ser virkelig grafen til $f(x)$ slik ut?} \textit{Var dette resultatet som forventet?} \textit{Kan dette skyldes numeriske feil?} 
\end{itemize}

\noindent
Husk at de to obligatoriske oppgavene i MAT-INF 1100 begge m氓 best氓s for 氓 kunne g氓 opp til endelig eksamen i emnet.

\paragraph{S酶knad om utsettelse av innleveringsfrist}Hvis man blir syk eller av andre grunner trenger 氓 s酶ke om utsettelse av innleveringsfristen m氓 du ta kontakt med studieadministrasjonen ved Matematisk institutt (7. et. Niels Henrik Abels hus, e-post: studieinfo@math.uio.no) i god tid f酶r innleveringsfristen. 
For fullstendige retningslinjer for innlevering av obligatoriske oppgaver se
\begin{center}
uio.no/studier/admin/obligatoriske-aktiviteter/mn-math-oblig.html
\end{center}

\section{\LaTeX}

Vi oppmuntrer til bruk av \LaTeX. Dette er et slags programmeringsspr氓k for 氓 skrive naturfaglige dokumenter som det kan v忙re greit 氓 l忙re seg f酶rst som sist.
Hjelp til \LaTeX{} kan du f氓 omtrent overalt, fra utallige forum p氓 verdensveven til gruppel忙rere over hele MatNat-fakultetet. 
Dersom du trenger hjelp til 氓 komme i gang vil du finne \LaTeX-versjonen av denne teksten samme sted som du finner PDF-filen. Vi vil ogs氓 legge ut en \LaTeX-mal med eksempler p氓 hvordan inkludere kode og figurer som du kan bruke til 氓 gj酶re obligene i MAT-INF1100 om du vil. 

\newpage

\begin{oppgaver}
\oppgave Ved hjelp av en GPS har vi m氓lt farten $v$ til et objekt som beveger seg. M氓lingene er gjort ved $N+1$ tidspunkter $(t_i)_{i=0}^N$ slik at resultatet er en f酶lge av tall-par $(t_i,v_i)_{i=0}^N$ der $v_i$ angir farten ved tidspunktet $t_i$.

\begin{deloppgaver}
\oppgave Gi en algoritme for 氓 beregne en tiln忙rming til objektets aksellerasjon $a(t)=v'(t)$ ut fra de beregnede verdiene $(t_i,v_i)$ av farten.

\oppgave Gi en algoritme for 氓 beregne en tiln忙rming til objektets avstand $s(t)$ fra startpunktet ut fra de beregnede verdiene n氓r $v(t)=s'(t)$ og $s(t_0)=0$.

\oppgave Fila 
\begin{center}
%\small
\url{http://www.uio.no/studier/emner/matnat/math/MAT-INF1100/h18/obliger/running.txt} 
\end{center}
er en logfil fra en l酶petur, der vi p氓 hver linje finner kommaseparerte tid/fart-verdier. Du har l忙rt at du kan lese inn verdiene fra denne fila i to
vektorer {\tt t} og {\tt v} ved hjelp av f酶lgende kode:
\begin{lstlisting}
t = []
v = []
infile = open('running.txt','r')
for line in infile:
    tnext, vnext = line.strip().split(',')
    t.append(float(tnext))
    v.append(float(vnext))
infile.close()
\end{lstlisting}
Last ned fila {\tt running.txt} og kj酶r denne koden, og bruk algoritmene fra a) og b) til 氓 lage to plott: Ett der du plotter objektets akselerasjon mot tid, og ett der du plotter
objektets avstand fra startpunktet mot tid.
\end{deloppgaver}





\oppgave 
Vi betrakter en pendel som ved tiden $t$ har et utslag p氓 $\theta(t)$ radianer (se Figur \ref{fig:pendel}). Ved hjelp av Newtons 2.~lov og litt vektorregning, kan vi utlede differensialligningen
\begin{equation}\label{eq:pendel}
\frac{d^2\theta}{dt^2}(t) = -\frac{g}{L}\sin(\theta(t))
\end{equation}
der $g \approx 9.81 \frac{m}{s^2}$ er gravitasjonsaksellerasjonen og $L$ er lengden til pendelen. Pendelligningen \eqref{eq:pendel} er et eksempel p氓 en enkel differensialligning som ikke har en l酶sning som uttrykkes ved enkel formel. I denne oppgaven skal vi se p氓 to m氓ter 氓 tiln忙rme l酶sningen p氓: \emph{Linearisering} og \emph{numerisk simulering}.
\begin{figure}
\centering
\begin{tikzpicture}[scale=2.5]
\draw[ultra thick] (0,0) -- (1,-1);
\draw[decorate,decoration={brace,amplitude=10pt},xshift=3pt,yshift=3pt] (0,0) -- (1.02,-1.02) node[midway,xshift=14pt,yshift=14pt]{$L$};
\draw[dashed] (0,0) -- (0,-1.3);
\draw[ball color=red] (1,-1) circle (3pt);
\begin{scope}
	\clip (0,0) -- (1,-1) -- (0,-1);
	\draw (0,0) circle (0.8);
\end{scope}
\draw (0.35,-0.82) node{$\theta$};
\end{tikzpicture}
\caption{En pendel best氓r av en stang med et tung objekt i den ene enden, og som dreier om et punkt i den andre enden. Vi antar at stangen er helt stiv, er tiln忙rmet massel酶s, og har lengde $L$.}
\label{fig:pendel}
\end{figure}

\begin{deloppgaver}
\oppgave
Hastigheten til pendelvekten er gitt ved $v(t) = L\frac{d\theta}{dt}(t)$. Forklar hvordan vi kan gj酶re om den \emph{annenordens} differensialligningen \eqref{eq:pendel} til \emph{systemet av f酶rsteordens} differensialligninger
\begin{equation}\label{eq:pendelsys}
\begin{split}
\frac{dv}{dt} &= -g\sin(\theta), \\[2pt]
\frac{d\theta}{dt} &= \frac{1}{L}v.
\end{split}
\end{equation}
Motsatt, hvordan finner du l酶sningen av \eqref{eq:pendel} om du kjenner l酶sningen av  \eqref{eq:pendelsys}?

\oppgave
Forklar hvorfor det \emph{line忙re} systemet
\begin{equation}\label{eq:pendelsyslin}
\begin{split}
\frac{dv}{dt} &= -g\theta ,\\[2pt]
\frac{d\theta}{dt} &= \frac{v}{L}
\end{split}
\end{equation}
gir en god tiln忙rming til l酶sningen av det \emph{ikkeline忙re} systemet \eqref{eq:pendelsys} n氓r utslaget til pendelen er lite. Finn l酶sningen av \eqref{eq:pendelsyslin} gitt startbetingelser $\theta(0)=\theta_0$, $v(0)=v_0$.

\noindent
\textit{Hint: For 氓 l酶se \eqref{eq:pendelsyslin}, skriv om ligningen til en annenordens differensialligning ved 氓 derivere ligningen for $\theta$.}

\oppgave
For et lite tall $h>0$ lar vi $t_k = hk$ for $k=0,1,\dots,N$, og vi lar $\theta_k, v_k$ v忙re tiln忙rminger til $\theta(t_k), v(t_k)$. Eulers metode for \eqref{eq:pendelsyslin} er da
\begin{equation}\label{eq:eulerpendellin}
\begin{split}
v_{k+1} &= v_k - gh\theta_k, \\
\theta_{k+1} &= \theta_k + h\frac{v_k}{L},
\end{split}
\end{equation}
for $k=0$, 1, \dots, $N-1$.
Skriv en Python-funksjon \verb|lin_pendel_euler| som med funksjonskallet 
\begin{verbatim}
v, theta = lin_pendel_euler(v0, theta0, g, L, N, h);
\end{verbatim}
bruker Eulers metode for 氓 beregne $v$ og $\theta$ gitt startbetingelser $v_0=\texttt{v0}, \theta_0=\texttt{theta0}$. Returverdiene \texttt{v} og \texttt{theta} skal v忙re Python-lister av lengde $N+1$ og inneholde verdiene $\texttt{v[0]}=v_0, \texttt{v[1]}=v_1, \dots, \texttt{v[N]}=v_N$, og tilsvarende for \texttt{theta}.

\oppgave
Plott den numeriske l酶sningen for vinkelutslaget som en funksjon av tid sammen med den eksakte l酶sningen du fant i b). Kommenter forskjellene. Bruk verdiene $g=9.81$, $L=1$, $v_0=0$, $\theta_0=\frac{\pi}{2}$ og $h = T/N$, der $T=4$. Bruk f酶rst $N=2^{5}$ antall steglengder og deretter $N=2^{10}$.

\oppgave
For 氓 beregne feilen i Eulers metode definerer vi 
\[
\epsilon(h) = |\theta_N-\theta(t_N)|,
\]
der $\theta_N$ er beregnet med Eulers metode \eqref{eq:eulerpendellin} med steglengde $h$, og $\theta(t)$ er l酶sningen du fant i b). (Vi ser bort fra feilen i $v$ her.) Beregn feilen $\epsilon(h)$ med verdiene $h=T/N$, der $N=2^{4},2^{5},\dots,2^{10}$. Bruk ellers de samme verdiene som i oppgave d). Oppgi feilene $\epsilon(h)$ i en tabell.

\oppgave
Vi sier at en numerisk metode er av orden $p$ dersom det finnes en konstant $C>0$ slik at $\epsilon(h) \leq Ch^p$ for sm氓 verdier av $h$. Vi kan finne tallet $p$ \emph{eksperimentelt} som f酶lger: Anta at $\epsilon(h)=Ch^p$ for verdier av $C$ og $p$ (som vi ikke kjenner). Velg to forskjellige verdier $h=h_1,h_2$. Da er $\epsilon(h_1)/\epsilon(h_2) = (Ch_1^p)/(Ch_2^p) = (h_1/h_2)^p$. L酶ser vi for $p$, f氓r vi
\[
p = \frac{\log\big(\epsilon(h_1)/\epsilon(h_2)\big)}{\log(h_1/h_2)}.
\]
Bruk tallene du fant i oppgave e) til 氓 estimere konvergensraten $p$ for verdiene $h_1=T/N$ og $h_2=h_1/2$, der $N=2^{4},2^5\dots,2^{9}$. Oppgi estimatene for $p$ i en tabell. Hva tror du konvergensraten til Eulers metode er?

\oppgave
Dersom utslaget til pendelen er stort, vil l酶sningen av det lineariserte systemet \eqref{eq:pendelsyslin} raskt avvike fra l酶sningen av den opprinnelige ligningen \eqref{eq:pendelsys}. Skriv en Python-funksjon \verb|pendel_euler| som med funksjonskallet 
\begin{verbatim}
v, theta = pendel_euler(v0, theta0, g, L, N, h);
\end{verbatim}
bruker Eulers metode for 氓 l酶se \eqref{eq:pendelsys}. Plott $\theta$ som en funksjon av tid fram til tiden $T=4$. Bruk en steglengde $h=T/N$, der $N=2^{10}$, og bruk ellers de samme verdiene som i oppgave d). Inkluder de to l酶sningene fra b) og c) i plottet. Kommenter forskjellene.

\oppgave
\emph{(Hvis du har tid)} \\
Systemet \eqref{eq:pendelsys} er et \emph{hamiltonsk system}, som kjennetegnes bl.a.~ved at det finnes en funksjon $H(v,\theta)$, \emph{hamiltonfunksjonen}, som er slik at $H(v(t),\theta(t))$ er konstant over tid. I dette tilfellet er $H(v,\theta)=-g\cos(\theta)+v^2/(2L)$, som kan tolkes som total (potensiell + kinetisk) energi for pendelen.
\begin{enumerate}
\oppgave Gjenta oppgave g), men med $T=50$. (N氓 trenger du ikke plotte l酶sningene av det line忙re systemet.) Hvordan tolker du plottet? Er energien til pendelen bevart i den numeriske l酶sningen? Hva skjer n氓r du velger en st酶rre $N$?
\oppgave Eulers metode har en tendens til 氓 tilf酶re energi til hamiltonske systemer. Skriv en Python-funksjon \verb|pendel_symp_euler| som implementerer den \emph{symplektiske Euler-metoden}\footnote{Symplektiske metoder er laget for 氓 bevare hamiltonfunksjonen langt mer n酶yaktig enn andre metoder. Disse metodene er sv忙rt viktige for 氓 simulere hamiltonske systemer slik som planetbaner og partikler i partikkelakselleratorer.}
\begin{equation}\label{eq:sympeulerpendellin}
\begin{split}
v_{k+1} &= v_k - gh\sin(\theta_k) \\
\theta_{k+1} &= \theta_k + h\frac{v_{k+1}}{L}.
\end{split}
\end{equation}
Gjenta oppgave g) med denne metoden, men med $T=50$ og med $T=500$. Ser den numeriske l酶sningen riktigere ut med denne metoden?
\end{enumerate}
\end{deloppgaver}





\vspace{1cm}
\centerline{\em Lykke til!}

\end{oppgaver}




\end{document}




