<!--
HTML file automatically generated from DocOnce source
(https://github.com/doconce/doconce/)
doconce format html mydoc.do.txt
-->
<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<meta name="generator" content="DocOnce: https://github.com/doconce/doconce/" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" />
<title></title>
<style type="text/css">
/* blueish style */
body {
margin-top: 1.0em;
background-color: #ffffff;
font-family: Helvetica, Arial, FreeSans, san-serif;
color: #000000;
}
h1 { font-size: 1.8em; color: #1e36ce; }
h2 { font-size: 1.6em; color: #1e36ce; }
h3 { font-size: 1.4em; color: #1e36ce; }
h4 { font-size: 1.2em; color: #1e36ce; }
a { color: #1e36ce; text-decoration:none; }
tt { font-family: "Courier New", Courier; }
p { text-indent: 0px; }
hr { border: 0; width: 80%; border-bottom: 1px solid #aaa}
p.caption { width: 80%; font-style: normal; text-align: left; }
hr.figure { border: 0; width: 80%; border-bottom: 1px solid #aaa; }
/* pre style removed because it will interfer with pygments */
div.highlight {
border: 1px solid #cfcfcf;
border-radius: 2px;
line-height: 1.21429em;
}
div.cell {
width: 100%;
padding: 5px 5px 5px 0;
margin: 0;
outline: none;
}
div.input {
page-break-inside: avoid;
box-orient: horizontal;
box-align: stretch;
display: flex;
flex-direction: row;
align-items: stretch;
}
div.inner_cell {
box-orient: vertical;
box-align: stretch;
display: flex;
flex-direction: column;
align-items: stretch;
box-flex: 1;
flex: 1;
}
div.input_area {
border: 1px solid #cfcfcf;
border-radius: 4px;
background: #f7f7f7;
line-height: 1.21429em;
}
div.input_area > div.highlight {
margin: .4em;
border: none;
padding: 0;
background-color: transparent;
}
div.output_wrapper {
position: relative;
box-orient: vertical;
box-align: stretch;
display: flex;
flex-direction: column;
align-items: stretch;
}
.output {
box-orient: vertical;
box-align: stretch;
display: flex;
flex-direction: column;
align-items: stretch;
}
div.output_area {
padding: 0;
page-break-inside: avoid;
box-orient: horizontal;
box-align: stretch;
display: flex;
flex-direction: row;
align-items: stretch;
}
div.output_subarea {
padding: .4em .4em 0 .4em;
box-flex: 1;
flex: 1;
}
div.output_text {
text-align: left;
color: #000;
line-height: 1.21429em;
}
div { text-align: justify; text-justify: inter-word; }
.tab {
padding-left: 1.5em;
}
div.toc p,a {
line-height: 1.3;
margin-top: 1.1;
margin-bottom: 1.1;
}
</style>
</head>
<!-- tocinfo
{'highest level': 1,
'sections': [("Solve the world's simplest differential equation",
1,
None,
'solve-the-world-s-simplest-differential-equation'),
('Mathematical problem', 2, None, 'mathematical-problem'),
('Numerical solution method',
2,
None,
'numerical-solution-method'),
('Implementation', 2, None, 'implementation'),
('Numerical experiments', 2, None, 'numerical-experiments')]}
end of tocinfo -->
<body>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({
TeX: {
equationNumbers: { autoNumber: "AMS" },
extensions: ["AMSmath.js", "AMSsymbols.js", "autobold.js", "color.js"]
}
});
</script>
<script type="text/javascript" async
src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.1/MathJax.js?config=TeX-AMS-MML_HTMLorMML">
</script>
<!-- ------------------- main content ---------------------- -->
<center>
<h4>Jun 22, 2021</h4>
</center> <!-- date -->
<br>
<h1 id="solve-the-world-s-simplest-differential-equation">Solve the world's simplest differential equation </h1>
<h2 id="mathematical-problem">Mathematical problem </h2>
<p>This exercise addresses the differential equation problem</p>
$$
\begin{align}
u'(t) &= -au(t), \quad t \in (0,T], \label{ode}\\
u(0) &= I, \label{initial:value}
\end{align}
$$
<p>where \( a \), \( I \), and \( T \) are prescribed constant parameters, and \( u(t) \) is
the unknown function to be estimated. This mathematical model
is relevant for physical phenomena featuring exponential decay
in time.
</p>
<h2 id="numerical-solution-method">Numerical solution method </h2>
<p>Derive the \( \theta \)-rule scheme for solving \eqref{ode} numerically
with time step \( \Delta t \):
</p>
$$
u^{n+1} = \frac{1 - (1-\theta) a\Delta t}{1 + \theta a\Delta t}u^n,
$$
<p>Here, \( n=0,1,\ldots,N-1 \).</p>
<b>Hint.</b>
<p>Set up the Forward Euler, Backward Euler, and the Crank-Nicolson (or
Midpoint) schemes first. Then generalize to the \( \theta \)-rule above.
</p>
<h2 id="implementation">Implementation </h2>
<p>The numerical method is implemented in a Python function
<code>solver</code> (found in the <a href="https://github.com/doconce/INF5620/tree/gh-pages/src/decay/experiments/decay_mod.py" target="_self"><tt>decay_mod</tt></a> module):
</p>
<!-- code=python (!bc pycod) typeset with pygments style "default" -->
<div class="cell border-box-sizing code_cell rendered">
<div class="input">
<div class="inner_cell">
<div class="input_area">
<div class="highlight" style="background: #f8f8f8">
<pre style="line-height: 125%;"><span style="color: #008000; font-weight: bold">from</span> <span style="color: #0000FF; font-weight: bold">numpy</span> <span style="color: #008000; font-weight: bold">import</span> linspace, zeros
<span style="color: #008000; font-weight: bold">def</span> <span style="color: #0000FF">solver</span>(I, a, T, dt, theta):
<span style="color: #BA2121; font-style: italic">"""Solve u'=-a*u, u(0)=I, for t in (0,T] with steps of dt."""</span>
dt <span style="color: #666666">=</span> <span style="color: #008000">float</span>(dt) <span style="color: #408080; font-style: italic"># avoid integer division</span>
N <span style="color: #666666">=</span> <span style="color: #008000">int</span>(<span style="color: #008000">round</span>(T<span style="color: #666666">/</span>dt)) <span style="color: #408080; font-style: italic"># no of time intervals</span>
T <span style="color: #666666">=</span> N<span style="color: #666666">*</span>dt <span style="color: #408080; font-style: italic"># adjust T to fit time step dt</span>
u <span style="color: #666666">=</span> zeros(N<span style="color: #666666">+1</span>) <span style="color: #408080; font-style: italic"># array of u[n] values</span>
t <span style="color: #666666">=</span> linspace(<span style="color: #666666">0</span>, T, N<span style="color: #666666">+1</span>) <span style="color: #408080; font-style: italic"># time mesh</span>
u[<span style="color: #666666">0</span>] <span style="color: #666666">=</span> I <span style="color: #408080; font-style: italic"># assign initial condition</span>
<span style="color: #008000; font-weight: bold">for</span> n <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(<span style="color: #666666">0</span>, N): <span style="color: #408080; font-style: italic"># n=0,1,...,N-1</span>
u[n<span style="color: #666666">+1</span>] <span style="color: #666666">=</span> (<span style="color: #666666">1</span> <span style="color: #666666">-</span> (<span style="color: #666666">1-</span>theta)<span style="color: #666666">*</span>a<span style="color: #666666">*</span>dt)<span style="color: #666666">/</span>(<span style="color: #666666">1</span> <span style="color: #666666">+</span> theta<span style="color: #666666">*</span>dt<span style="color: #666666">*</span>a)<span style="color: #666666">*</span>u[n]
<span style="color: #008000; font-weight: bold">return</span> u, t
</pre>
</div>
</div>
</div>
</div>
<div class="output_wrapper">
<div class="output">
<div class="output_area">
<div class="output_subarea output_stream output_stdout output_text">
</div>
</div>
</div>
</div>
</div>
<h2 id="numerical-experiments">Numerical experiments </h2>
<p>Fix the values of where \( I \), \( a \), and \( T \).
Then vary \( \Delta t \) for \( \theta=0,1/2,1 \).
Illustrate that if \( \Delta t \) is not sufficiently small,
\( \theta=0 \) and \( \theta=1/2 \) can give non-physical solutions
(more precisely, oscillating solutions).
</p>
<p>Perform experiments and determine empirically the convergence
rate for \( \theta=0,1/2,1 \).
</p>
<!-- ------------------- end of main content --------------- -->
</body>
</html>