\input{../preamble}

\title{POTETEN}

\begin{document}
\pagenumbering{gobble}

%primtall

\maketitle


Emnepakken TMA4101/4106/4111/4121
har som grunnfilosofi at innholdet skal motiveres fra anvendelser.
I hvert fag er det mulig å velge et langsgående prosjekt som obligatorisk arbeid.
I hvert emne kan du velge \textbf{standardprosjektet},
som vi omhandle varmeledning.
Poenget med prosjektene er å trene på praktisk modellering,
samt drive gjennom at en modells 
kvalitet måles ved å sjekke om den gir korrekte prediksjoner ved empirisk eksperiment,
All vitenskap er tuftet på dette. 


\section*{4101}


Dersom du har en kokt elgtunge med temperatur $T(t)$ som står og kjøles ned 
(eller varmes opp) i omgivelser med temperatur $T_k$,
er Newtons avkjølingslov
\[
\dot{T}(t)=\alpha \left(T(t)-T_K\right) \qquad T(0)=T_0
\]
en grov modell for avkjølingen (eller oppvarmingen).
Proporsjonalitetskonstanten $\alpha$ inneholder informasjon om 
varmekapasiteten til elgtungen og hvor fort
varmeflyten går mellom elgtungen og omgivelsene.
Løsningen til differensiallikningen er
\[
T(t) = T_K + (T_0-T_K)e^{- \alpha t}.
\]
I dette prosjektet skal vi teste modellen empirisk. 
Velg en gjenstand, 
for eksempel et glass vann eller en potet eller en kokt elgtunge eller noe helt annet, 
la den avkjøles mens du måler temperaturen,
og sammenlikne med modellen. 
I praksis blir jobben å finne korrekt $\alpha$ for ditt oppsett. 
Jeg gjorde dette med en rektangulær skål med vann,
beregnet $\alpha$ fra målingene,
og fikk figuren under. 
Som du ser er det antagelig noe fysikk her som Newtons avkjølingslov ikke tar høyde for,
for eksempel fordampning.

\begin{center}
\includegraphics[scale=.7]{../../figurer/vannet}\\
Blå: Målt temperatur \hspace{10mm}
Oransje: Newtons avkjølingslov
\end{center}

\clearpage

\section*{4106}








Når man skal løse sin første partielle differensiallikning numerisk,
er varmelikningen grei å begynne med, 
for dette er den enkleste å få til.
Den deriverte til en funksjon $f$ er
\[
f'(x)=\lim_{h \to 0} \frac{f(x+h)-f(x)}{h}
\]
dersom grensen eksisterer.
Uttrykket 
\[
\frac{f(x+h)-f(x)}{h}
\]
er stigningstallet til sekanten til $f$ mellom punktene $x$ og $x+h$. 
For små $h$ er denne sekanten en grei tilnærming til stigningstallet $f(x)$:
\[
f'(x)\approx \frac{f(x+h)-f(x)}{h}.
\]
På papiret er det slik at jo mindre $h$, desto bedre tilnærming. 
La $f(x)=e^x$, og la $h=0.1$. Vi beregner
\[
f'(1.5)\approx \frac{e^{1.6}-e^{1.5}}{0.1}=4.7134. 
\]
Merk at 
\[
f'(1.5)=e^{1.5}=4.4817,
\]
så denne tilnærmingen bommer med rundt $2 \cdot 10^{-1}$.

\begin{oppgave}{1}
Gjenta eksperimentet med  $h=0.01$, $h=0.001$ osv,
og se hva som skjer. 
Hvor liten $h$ kan du ta før det går åt skogen?
\end{oppgave}


Feilen i forrige  er tydelig proporsjonal med $h$ - 
deler du $h$ på 10, deler du feilen på 10. 
En god illustrasjon av lineær feil. 
Men dette gjelder bare inntil et visst punkt, 
blir $h$ liten nok blir presisjonen dårligere igjen. 
Dette er fordi datamaskinen din kun regner med 16 desimaler. 
Skal vi få til noe bedre, 
må vi bruke en annen formel.


Det går an å utlede formelen
\[
f'(x)\approx \frac{f(x+h)-f(x)}{h}.
\]
ved å bruke taylorutviklingen til $f$ i punktet $x$:
\[
f(x+h)=f(x)+f'(x)h+\frac{f''(x)}{2}h^2+\cdots.
\]
Stigningen til sekanten kan skrives
\[
\frac{f(x+h)-f(x)}{h}=f'(x)+\frac{f''(x)}{2}h+\cdots
\]
og vi ser at dette  stigninstallet består av den eksakte verdien
for $f'(x)$ pluss resten av taylorrekken til $f$
\[
\frac{f''(x)}{2}h+\frac{f'''(x)}{6}h^2+\cdots
\]
Denne halen forteller oss noe om feilen.
Dersom  $h$ er liten nok, vil $h$ være mye større enn $h^2$, og vi skriver
\[
\frac{f''(x)}{2}h+\frac{f'''(x)}{6}h^2+\cdots=O(h)
\]
for å signalisere at feilen er proporsjonal med $h$.


\begin{oppgave}{2}
Gjenta eksperimentet i forrige oppgave, 
men nå med formelen
\[
\frac{f(x+h)-f(x-h)}{2h}=f'(x)
\] 
og se hva som skjer. 
Hvor liten $h$ kan du nå ta før det går åt skogen?
Bruk taylorrekker til å forklare oppførselen.
\end{oppgave}


Legg merke til hvordan feilen i forrige eksempel er nærmest perfekt kvadratisk - 
vi får to nye desimaler hver gang vi deler $h$ på 10. Feilen deles altså på 100 når $h$ deles på 10.


\begin{oppgave}{3}
Hvis du virkelig vil slå på stortrommen, kan du bruke formelen
\[
f'(x)= \frac{f(x-2h)-8f(x-h)+8f(x+h)-f(x+2h)}{12h}
\]
Prøv.
\end{oppgave}


Dersom du trenger en tilnærming for $f''(x)$, 
kan du bruke den andre ordens sentraldifferansen
\[
\frac{f(x+h)-2f(x)+f(x-h)}{h^2}=f''(x)+\frac{f''''(x)}{12}h^2+\cdots
\]
Det finnes bedre og mer kompliserte formler, 
men vi skal ikke plages med dem.
Vi skal nå lage en numerisk metode for varmeligningen
\[
\dot u(x,t)=u''(x,t)
\]
med randkrav
\[
u(0,t)=u(1,t)=0
\]
og initialkrav
\[
u(x,0)=f(x).
\]
Det er relativt lett å konstruere liknende metoder for andre likninger
men varmelikningen er som sagt den enkleste å begynne med.
Når vi skal løse en partiell differensiallikning, 
må vi holde styr på to gitre - et i $x$-retningen, og et i $t$-retningen.

\begin{tcolorbox}
Vi skal finne approksimasjoner til løsningen $u$ i alle gitterpunktene.
\end{tcolorbox}

Vi gitrer opp intervallet $[0,1]$ på $x$-aksen med gitteravstanden $h$, 
og nummererer punktene slik at $x_0=0$ og $x_n=1$.
Den positive $t$-aksen gitrer vi opp med gitteravstanden $k$, 
og nummerer slik at $t_0=0$.
\begin{center}
\hspace*{-3em}
    \begin{tikzpicture}[thick]
%      \draw[->] (-1.5,0) -- (1.5,0) node[right]{$x$};
          \draw[->] (-3,-1) -- (4,-1) node[right]{$x$};
          \draw[->] (-5*.6,-1) -- (-5*.6,3) node[right]{$t$};
           \draw[-] (5*.6,-1) -- (5*.6,3) ;

%          \draw[->] (-2.5,-1) -- (-2.5,1.3) node[above]{$y$};
%          \draw[line width = 1mm] (-2.35,.05) -- (-2.35,-.05) node[below]{$a$};
          \draw[line width = 1mm] (5*.6,-.95) -- (5*.6,-1.05) node[below]{$x=1$};
           \draw[line width = 1mm] (4*.6,-.95) -- (4*.6,-1.05) ;
           \draw[line width = 1mm] (3*.6,-.95) -- (3*.6,-1.05) ;
           \draw[line width = 1mm] (2*.6,-.95) -- (2*.6,-1.05) ;
           \draw[line width = 1.4mm] (1*.6,-.93) -- (1*.6,-1.07) node[below]{$x_i=ih$};
           \draw[line width = 1.4mm] (-5*.586,1) -- (-5*.614,1) node[left]{$t_j=jk$};
           \draw[line width = 1.4mm] (1*.586,1) -- (1*.614,1) node[above]{$(x_i,t_j)=(ih,jk)$};
           \draw[line width = 1mm] (0,-.95) -- (0,-1.05) ;
           \draw[line width = 1mm] (-1*.6,-.95) -- (-1*.6,-1.05) ;
           \draw[line width = 1mm] (-2*.6,-.95) -- (-2*.6,-1.05) ;
           \draw[line width = 1mm] (-3*.6,-.95) -- (-3*.6,-1.05) ;
           \draw[line width = 1mm] (-4*.6,-.95) -- (-4*.6,-1.05) ;        
           \draw[line width = 1mm] (-5*.6,-.95) -- (-5*.6,-1.05) ;
          \draw[line width = 1mm] (5*.6,.05) -- (5*.6,-.05) ;
           \draw[line width = 1mm] (4*.6,.05) -- (4*.6,-.05) ;
           \draw[line width = 1mm] (3*.6,.05) -- (3*.6,-.05) ;
           \draw[line width = 1mm] (2*.6,.05) -- (2*.6,-.05) ;
           \draw[line width = 1mm] (1*.6,.05) -- (1*.6,-.05) ;
           \draw[line width = 1mm] (0,.05) -- (0,-.05) ;
           \draw[line width = 1mm] (-1*.6,.05) -- (-1*.6,-.05) ;
           \draw[line width = 1mm] (-2*.6,.05) -- (-2*.6,-.05) ;
           \draw[line width = 1mm] (-3*.6,.05) -- (-3*.6,-.05) ;
           \draw[line width = 1mm] (-4*.6,.05) -- (-4*.6,-.05) ;        
           \draw[line width = 1mm] (-5*.6,.05) -- (-5*.6,-.05) ;

          \draw[line width = 1mm] (5*.6,.95) -- (5*.6,1.05) ;
           \draw[line width = 1mm] (4*.6,.95) -- (4*.6,1.05) ;
           \draw[line width = 1mm] (3*.6,.95) -- (3*.6,1.05) ;
           \draw[line width = 1mm] (2*.6,.95) -- (2*.6,1.05) ;
           \draw[line width = 1.4mm] (1*.6,.93) -- (1*.6,1.07) ;
           \draw[line width = 1mm] (0,.95) -- (0,1.05) ;
           \draw[line width = 1mm] (-1*.6,.95) -- (-1*.6,1.05) ;
           \draw[line width = 1mm] (-2*.6,.95) -- (-2*.6,1.05) ;
           \draw[line width = 1mm] (-3*.6,.95) -- (-3*.6,1.05) ;
           \draw[line width = 1mm] (-4*.6,.95) -- (-4*.6,1.05) ;        
           \draw[line width = 1mm] (-5*.6,.95) -- (-5*.6,1.05) ;

          \draw[line width = 1mm] (5*.6,1.95) -- (5*.6,2.05) ;
           \draw[line width = 1mm] (4*.6,1.95) -- (4*.6,2.05) ;
           \draw[line width = 1mm] (3*.6,1.95) -- (3*.6,2.05) ;
           \draw[line width = 1mm] (2*.6,1.95) -- (2*.6,2.05) ;
           \draw[line width = 1mm] (1*.6,1.95) -- (1*.6,2.05) ;
           \draw[line width = 1mm] (0,1.95) -- (0,2.05) ;
           \draw[line width = 1mm] (-1*.6,1.95) -- (-1*.6,2.05) ;
           \draw[line width = 1mm] (-2*.6,1.95) -- (-2*.6,2.05) ;
           \draw[line width = 1mm] (-3*.6,1.95) -- (-3*.6,2.05) ;
           \draw[line width = 1mm] (-4*.6,1.95) -- (-4*.6,2.05) ;        
           \draw[line width = 1mm] (-5*.6,1.95) -- (-5*.6,2.05) ;


            \draw[line width = 1mm] (-5*.59,0) -- (-5*.61,0) ;
            \draw[line width = 1mm] (-5*.59,2) -- (-5*.61,2) ;
            \end{tikzpicture}
  \end{center}


Vi setter først opp en tilnærming for $u''(x,t)$, 
basert på den andre ordens endelige differanseformelen for $x$:
\[
u''(x,t) \approx \frac{u(x+h,t)-2u(x,t)+u(x-h,t)}{h^2}
\]
Vi kan nå tenke at vi erstatter 
$u(x,t)$ med $n+1$ envariable funksjoner $u_i(t)$, 
som beskriver temperaturendringen i hvert sitt punkt $x_i$ på stangen:
\[
u''(x_i,t) \approx \frac{u_{i+1}(t)-2u_{i}(t)+u_{i-1}(t)}{h^2},
\]
Nå setter vi inn dette uttrykket i varmelikningen: 
\[
\dot u(x_i,t)=u''(x_i,t) \approx \frac{u_{i+1}(t)-2u_{i}(t)+u_{i-1}(t)}{h^2},
\]
Siden $u_i$ kun er en funksjon av $t$, 
passer det å skrive alt dette som et koblet system av differensiallikninger: 
\[
\dot u_i(t) = \frac{u_{i+1}(t)-2u_{ij}(t)+u_{i-1}(t)}{h^2}
\]
Systemet har $n-1$ likninger, 
for $u_0=0$ og $u_n=0$ kjenner vi; 
disse er gitt av randbetingelsene.


Nå kan vi løse systemet med en ønsket numerisk metode. 
og vi skal ta i betraktning tre valg. 
Etter diskretiseringen i $t$, 
skriver vi approksimasjonen i punktet $(x_i,t_j)$ som $u_{ij}$:
\[
u(x_i,t_j)\approx u_{ij}
\]
De tre metodene vi skal bruke, 
er eksplisitt Euler, 
implisitt Euler, 
og trapesmetoden. 
De korresponderende skjemaene for varmelikningen blir: \\

\textbf{Eksplisitt}\hspace{20mm}
$
\displaystyle
\frac{u_{i,j+1}-u_{ij}}{k} =\frac{u_{i+1,j}-2u_{i,j}+u_{i-1,j}}{h^2}  \qquad \quad 
$

\textbf{Imsplisitt}\hspace{20mm}
$
\displaystyle
\frac{u_{i,j+1}-u_{ij}}{k} =\frac{u_{i+1,j+1}-2u_{i,j+1}+u_{i-1,j+1}}{h^2} 
$

\textbf{Crank-Nicolson}\hspace{10mm}
$
\displaystyle
\frac{u_{i,j+1}-u_{ij}}{k} = \frac{u_{i+1,j}-2u_{i,j}  + u_{i-1,j}}{2h^2} +\frac{u_{i+1,j+1}-2u_{i,j+1}+u_{i-1,j+1}}{2h^2}
$

Metodene har litt forskjellige egenskaper. 
Vi får se hva vi rekker av analyse.

\begin{oppgave}{4}
Løs varmelikningen med eksplisitt skjema. 
Prøv forskjellige kombinasjoner av $h$ og $k$.\\
(Metoden funker ikke med mindre $k/h^2\geq 1/2$.)
\end{oppgave}

\begin{oppgave}{5}
Kravet $k/h^2\geq 1/2$ er ganske restriktivt,
og det er derfor vi har implisitte metoder. 
De er litt vanskeligere å implementere, 
og litt mer beregningstunge,
men slipper unna slike hårreisende betingelser på $h$ og $k$.
Implementer implisitt.
\end{oppgave}


\begin{oppgave}{6}
Men hvorfor har vi Crank-Nicolson?
Implementer Crank-Nicolson, 
og sammenlikne de tre metodene for samme $h$ og $k$ med den 
analytiske løsningen til varmelikningen med for eksempel initialkrav $f(x)=\sin x$.
\end{oppgave}



\clearpage

\section*{4111}


Et av to store lys i tunnelen dette semesteret er varmelikningen\footnote{\url{https://en.wikipedia.org/wiki/Heat_equation}}
\[
\dot u = \alpha \Delta u
\]
der $\alpha$  kalles den termiske diffusiviteten.\footnote{\url{https://en.wikipedia.org/wiki/Thermal_diffusivity}}. 
I dette semesteret skal vi prøve å estimere den termiske diffusiviteten til en typisk potet.
Du er nødt til å skaffe et par poteter og et steketermometer,
koke eller steke et par poteter og sammenlikne med analytisk eller numerisk løsning av varmelikningen. 
Du må sette opp en modell og gjøre forenklende antakelser. 
Hvis du klarer å finne termisk diffusivitet, 
kan prøve å gå videre og finne varmekapasitet og termisk konduktivitet. 

Du kan også bruke noe annet enn en potet,
for eksempel en lammestek.
Hvis du har forslag til andre prosjekter, 
er det bare å foreslå,
send meg en epost.
For eksempel kan Lorang og kompani fortsette med potensiostaten sin,
og Jenny kan herje videre med den analoge datamaskinen. \\

\includegraphics[scale=.465]{../../figurer/sprengt-potet-1}
\begin{center}
En gang sprengte en potet i ovnen.
\end{center}

\clearpage

Jeg har nå begynt på dette prosjektet selv. 
De termometrene jeg hadde til rådighet i Vallegata 8B (her i Norges hippeste strøk St. Hanshaugen, 
naboene og foreldrene i barnehagen og så videre er berømtheter alle sammen) 
var ikke gode nok til å måle temperaturen i en potet. 
Steketermometeret var for stort (varme fra omgivelsene lekket inn og ødela alt),
og utetermometeret vårt likte ikke så godt å måle temperatur inni en potet. 
Så om du skal få dette til med en potet må du ha et bedre termometer enn jeg har.
Jeg ønsker meg et bedre termometer til jul. 

Men du kan som sagt estimere den termiske diffusiviteten til akkurat hva du vil. 
Jeg landa til slutt på en rektangulær skål med vann,
som jeg lot avkjøle i vinduskarmen. 
Termometeret var sentrert sånn omtrent i midt i boksen,
og målingene (oransje) er plottet i samme figur som den analytiske løsningen (blå).
Som du ser er ikke matchen så veldig god, 
men jeg har en ide om hvorfor,
og skal skrive en utførlig beskrivelse av hva jeg har gjort så snart jeg får tid. 
Jeg har gjort mange idiotiske feil og lært veldig mye,
og er ikke spesielt imponert over mitt eget arbeid,
men ville nok godkjent mitt eget prosjekt hvis jeg hadde levert det inn. 


\begin{minipage}{90mm}
\includegraphics[scale=.08]{../../figurer/img_0571}
\end{minipage}
\begin{minipage}{80mm}
\includegraphics[scale=.4]{poteten-II-analytisk}
\end{minipage}

Jeg herja noe voldsomt med penn og papir og hosta til slutt opp denne analytiske formelen:

\begin{center}
\includegraphics[scale=.45]{../../figurer/varmeboks-analytisk}
\end{center}

Jeg vet at poteten ikke er en rektangulær boks, 
men i fysikk er det tross alt noe som heter \say{assume a spherical cow}:\\
\url{https://en.wikipedia.org/wiki/Spherical_cow}\\
Hvis du har målt dimensjonene på det du målte temperaturen på og vet posisjonen til termometeret,
kan du bruke formelen over til å regne ut hvordan varmelikningen påstår at temperaturutviklingen skal se ut for en fritt valgt $\mathbf x$.

Jeg har implementert dette i \say{poteten-II-analytisk.py} her: \\
\url{https://folk.ntnu.no/mortano/prosjekt/poteten/}\\
og i rust her: \\
\url{https://folk.ntnu.no/mortano/rust/poteten/}\\

Pythonkoden er skrevet med for-løkker og kjører ganske treigt. 
Hvis noen gidder å skrive den om slik at den går fortere er det fritt frem. 
Rustkoden er også basert på for-løkker, 
men i motsetning til i python går dette unna, 
så du kan regne ut mange flere ledd i fourierrekken. 
Rustkoden produserer en .csv-fil som du kan flytte til en egnet mappe og plotte på et eller annet vis. 




\end{document}
