aboutsummaryrefslogtreecommitdiff
path: root/src/Lecture7/slides/.ipynb_checkpoints
diff options
context:
space:
mode:
authorSebastiano Tronto <sebastiano.tronto@gmail.com>2021-05-25 17:10:49 +0200
committerSebastiano Tronto <sebastiano.tronto@gmail.com>2021-05-25 17:10:49 +0200
commitd6c61d988bfa4255baf9cdae42db59ebee38363f (patch)
tree118ff3c2424e735149c145524965a4a337e50beb /src/Lecture7/slides/.ipynb_checkpoints
parent46eef66b1e1571c77dc828d7e950b129b4c8bfd0 (diff)
downloadmathsoftware-d6c61d988bfa4255baf9cdae42db59ebee38363f.tar.gz
mathsoftware-d6c61d988bfa4255baf9cdae42db59ebee38363f.zip
Added files
Diffstat (limited to 'src/Lecture7/slides/.ipynb_checkpoints')
-rw-r--r--src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-checkpoint.ipynb394
-rw-r--r--src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-notebook-checkpoint.ipynb412
-rw-r--r--src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-checkpoint.ipynb165
-rw-r--r--src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-notebook-checkpoint.ipynb298
4 files changed, 1269 insertions, 0 deletions
diff --git a/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-checkpoint.ipynb b/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-checkpoint.ipynb
new file mode 100644
index 0000000..ab7dd49
--- /dev/null
+++ b/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-checkpoint.ipynb
@@ -0,0 +1,394 @@
1{
2 "cells": [
3 {
4 "cell_type": "markdown",
5 "metadata": {},
6 "source": [
7 "# Nested loops\n",
8 "\n",
9 "The following two functions compute sum and product of matrices, respectively.\n",
10 "\n",
11 "By counting the nested loops it is easy to see that `add()` is $O(n^2)$ while `prod()` is $O(n^3)$."
12 ]
13 },
14 {
15 "cell_type": "code",
16 "execution_count": 37,
17 "metadata": {},
18 "outputs": [
19 {
20 "name": "stdout",
21 "output_type": "stream",
22 "text": [
23 "Time for add: 0.00012074300000008975\n",
24 "Time for prod: 0.00036587199999971176\n"
25 ]
26 }
27 ],
28 "source": [
29 "from random import randint\n",
30 "import time\n",
31 "\n",
32 "def add(A, B):\n",
33 " S = [[0] * len(A) for i in range(len(A))]\n",
34 " for i in range(len(A)):\n",
35 " for j in range(len(A)):\n",
36 " S[i][j] = A[i][j] + B[i][j]\n",
37 " return S\n",
38 "\n",
39 "def prod(A, B):\n",
40 " S = [[0] * len(A) for i in range(len(A))]\n",
41 " for i in range(len(A)):\n",
42 " for j in range(len(A)):\n",
43 " for k in range(len(A)):\n",
44 " S[i][j] = S[i][j] + A[i][k] * B[k][j]\n",
45 " return S\n",
46 "\n",
47 "N = 10\n",
48 "A = [ [randint(0,100) for i in range(N)] for j in range(N) ]\n",
49 "B = [ [randint(0,100) for i in range(N)] for j in range(N) ]\n",
50 "\n",
51 "t0 = time.process_time()\n",
52 "add(A,B)\n",
53 "t1 = time.process_time()\n",
54 "prod(A,B)\n",
55 "t2 = time.process_time()\n",
56 "\n",
57 "print(\"Time for add: \", t1-t0)\n",
58 "print(\"Time for prod:\", t2-t1)"
59 ]
60 },
61 {
62 "cell_type": "markdown",
63 "metadata": {},
64 "source": [
65 "# Sorting a list, slow version\n",
66 "\n",
67 "The following code implements a slow version of the so-called *insertion sort* alogithm\n",
68 "\n",
69 "Complexity: $O(n^2)$."
70 ]
71 },
72 {
73 "cell_type": "code",
74 "execution_count": 61,
75 "metadata": {},
76 "outputs": [
77 {
78 "name": "stdout",
79 "output_type": "stream",
80 "text": [
81 "Running time: 1.1012288430000012\n"
82 ]
83 }
84 ],
85 "source": [
86 "from random import randint\n",
87 "import time\n",
88 "\n",
89 "def correct_position(e, S):\n",
90 " for i in range(len(S)):\n",
91 " if S[i] > e:\n",
92 " return i\n",
93 " return len(S)\n",
94 "\n",
95 "def sort_list(L):\n",
96 " S = []\n",
97 " for e in L:\n",
98 " cp = correct_position(e, S)\n",
99 " S.insert(cp, e)\n",
100 " return S\n",
101 "\n",
102 "N = 10000\n",
103 "L = [randint(0,10**9) for i in range(N)]\n",
104 "\n",
105 "t0 = time.process_time()\n",
106 "sort_list(L)\n",
107 "t1 = time.process_time()\n",
108 "\n",
109 "print(\"Running time:\", t1-t0)"
110 ]
111 },
112 {
113 "cell_type": "markdown",
114 "metadata": {},
115 "source": [
116 "# Binary search\n",
117 "\n",
118 "The following code implements a binary search.\n",
119 "\n",
120 "Complexity: $O(\\log_2(n))$"
121 ]
122 },
123 {
124 "cell_type": "code",
125 "execution_count": 53,
126 "metadata": {},
127 "outputs": [
128 {
129 "name": "stdout",
130 "output_type": "stream",
131 "text": [
132 "The correct position of e = 658230309 in L is:\n",
133 "... 658211821 658224379 e 658234625 658246765 ...\n",
134 "\n",
135 "Time for sorting: 0.021211020999999164\n",
136 "Time for searching: 7.820199999741817e-05\n"
137 ]
138 }
139 ],
140 "source": [
141 "from random import randint\n",
142 "import time\n",
143 "\n",
144 "def binary_search(e, S, start, end):\n",
145 " if start == end:\n",
146 " return start\n",
147 " midpoint = (start+end) // 2\n",
148 " if e < S[midpoint]:\n",
149 " return binary_search(e, S, start, midpoint)\n",
150 " else:\n",
151 " return binary_search(e, S, midpoint+1, end)\n",
152 " \n",
153 "N = 100000\n",
154 "L = [randint(0,10**9) for i in range(N)]\n",
155 "e = randint(0,10**9)\n",
156 "\n",
157 "t0 = time.process_time()\n",
158 "L.sort() # Using Python's sort()\n",
159 "t1 = time.process_time()\n",
160 "i = binary_search(e, L, 0, len(L))\n",
161 "t2 = time.process_time()\n",
162 "print(\"The correct position of e =\", e, \"in L is:\")\n",
163 "print(\"...\", L[i-2], L[i-1], \"e\", L[i], L[i+1], \"...\")\n",
164 "print(\"\")\n",
165 "print(\"Time for sorting: \", t1-t0)\n",
166 "print(\"Time for searching:\", t2-t1)\n"
167 ]
168 },
169 {
170 "cell_type": "markdown",
171 "metadata": {},
172 "source": [
173 "# Sorting a list, fast version (with binary_search)\n",
174 "\n",
175 "The following code uses the function `binary_search()` above instead of `correct_position()` in our insertion sort algorithm.\n",
176 "\n",
177 "Complexity: $O(n\\log_2(n))$"
178 ]
179 },
180 {
181 "cell_type": "code",
182 "execution_count": 69,
183 "metadata": {},
184 "outputs": [
185 {
186 "name": "stdout",
187 "output_type": "stream",
188 "text": [
189 "Running time: 0.03710268399998995\n"
190 ]
191 }
192 ],
193 "source": [
194 "from random import randint\n",
195 "import time\n",
196 "\n",
197 "def binary_search(e, S, start, end):\n",
198 " if start == end:\n",
199 " return start\n",
200 " midpoint = (start+end) // 2\n",
201 " if e < S[midpoint]:\n",
202 " return binary_search(e, S, start, midpoint)\n",
203 " else:\n",
204 " return binary_search(e, S, midpoint+1, end)\n",
205 " \n",
206 "def sort_list(L):\n",
207 " S = []\n",
208 " for e in L:\n",
209 " cp = binary_search(e, S, 0, len(S)) # Changed here\n",
210 " S.insert(cp, e)\n",
211 " return S\n",
212 " \n",
213 "N = 10000\n",
214 "L = [randint(0,10**9) for i in range(N)]\n",
215 "\n",
216 "t0 = time.process_time()\n",
217 "sort_list(L)\n",
218 "t1 = time.process_time()\n",
219 "\n",
220 "print(\"Running time:\", t1-t0)"
221 ]
222 },
223 {
224 "cell_type": "markdown",
225 "metadata": {},
226 "source": [
227 "# Fast exponentiation\n",
228 "\n",
229 "The following cell contains two functions for computing $a^n$ ($n$ non-negative integer): a slow one that runs in $O(n)$ and a fast one that runs in $O(\\log_2(n))$. We compare these two also with Python's built-in operator `**`.\n",
230 "\n",
231 "Complexity: $O(n)$ for the slow algorithm, $O(\\log_2(n))$ for the other two."
232 ]
233 },
234 {
235 "cell_type": "code",
236 "execution_count": 30,
237 "metadata": {},
238 "outputs": [
239 {
240 "name": "stdout",
241 "output_type": "stream",
242 "text": [
243 "2.71828179834636\n",
244 "2.7182817863957984\n",
245 "2.7182817983473577\n",
246 "Time for slow_power(): 3.234879998000004\n",
247 "Time for fast_power(): 9.059099999575437e-05\n",
248 "Time for Python's **: 0.00010159500000384014\n"
249 ]
250 }
251 ],
252 "source": [
253 "import time\n",
254 "\n",
255 "def slow_power(a, n):\n",
256 " r = 1\n",
257 " for i in range(n):\n",
258 " r = r * a\n",
259 " return r\n",
260 "\n",
261 "def fast_power(a, n):\n",
262 " if n == 0:\n",
263 " return 1\n",
264 " if n%2 == 0:\n",
265 " return fast_power(a*a, n//2)\n",
266 " else:\n",
267 " return a * fast_power(a, n-1)\n",
268 "\n",
269 "a = 1.00000001\n",
270 "n = 100000000\n",
271 "\n",
272 "t0 = time.process_time()\n",
273 "print(slow_power(a, n))\n",
274 "t1 = time.process_time()\n",
275 "print(fast_power(a, n))\n",
276 "t2 = time.process_time()\n",
277 "print(a**n)\n",
278 "t3 = time.process_time()\n",
279 "\n",
280 "print(\"Time for slow_power():\", t1-t0)\n",
281 "print(\"Time for fast_power():\", t2-t1)\n",
282 "print(\"Time for Python's **: \", t3-t2)"
283 ]
284 },
285 {
286 "cell_type": "markdown",
287 "metadata": {},
288 "source": [
289 "# Fast gcd\n",
290 "\n",
291 "Complexity: $O(\\log_2(n))$"
292 ]
293 },
294 {
295 "cell_type": "code",
296 "execution_count": 31,
297 "metadata": {},
298 "outputs": [
299 {
300 "name": "stdout",
301 "output_type": "stream",
302 "text": [
303 "126\n",
304 "Running time: 0.00017707599999994272\n"
305 ]
306 }
307 ],
308 "source": [
309 "import time\n",
310 "\n",
311 "def gcd(a, b):\n",
312 " if b == 0:\n",
313 " return a\n",
314 " else:\n",
315 " return gcd(b, a%b)\n",
316 "\n",
317 "t0 = time.process_time()\n",
318 "print(gcd(155275387236018, 572335397352432))\n",
319 "t1 = time.process_time()\n",
320 "\n",
321 "print(\"Running time:\", t1-t0)"
322 ]
323 },
324 {
325 "cell_type": "markdown",
326 "metadata": {},
327 "source": [
328 "# Fibonacci numbers\n",
329 "\n",
330 "In the following cell there are two functions that compute the $n$-th Fibonacci number. They are almost the same, but the second one memorizes the results in a list to avoid computing them multiple times, and it is much much faster.\n",
331 "\n",
332 "Complexity: $O\\left(\\left(\\frac{1+\\sqrt 5}{2}\\right)^n\\right)\\sim O(1.6^n)$ for the slow version, $O(n)$ for the fast version."
333 ]
334 },
335 {
336 "cell_type": "code",
337 "execution_count": null,
338 "metadata": {},
339 "outputs": [],
340 "source": [
341 "import time\n",
342 "\n",
343 "F_memorized = [-1] * (10**6)\n",
344 "\n",
345 "def F_slow(n):\n",
346 " if n <= 1:\n",
347 " return n\n",
348 " else:\n",
349 " return F_slow(n-1) + F_slow(n-2)\n",
350 " \n",
351 "def F_fast(n):\n",
352 " if F_memorized[n] == -1:\n",
353 " if n <= 1:\n",
354 " F_memorized[n] = n\n",
355 " else:\n",
356 " F_memorized[n] = F_fast(n-1) + F_fast(n-2)\n",
357 " \n",
358 " return F_memorized[n]\n",
359 "\n",
360 "n = 40\n",
361 "\n",
362 "t0 = time.process_time()\n",
363 "print(F_slow(n))\n",
364 "t1 = time.process_time()\n",
365 "print(F_fast(n))\n",
366 "t2 = time.process_time()\n",
367 "\n",
368 "print(\"Time for F_slow:\", t1-t0)\n",
369 "print(\"Time for F_fast:\", t2-t1)"
370 ]
371 }
372 ],
373 "metadata": {
374 "kernelspec": {
375 "display_name": "Python 3",
376 "language": "python",
377 "name": "python3"
378 },
379 "language_info": {
380 "codemirror_mode": {
381 "name": "ipython",
382 "version": 3
383 },
384 "file_extension": ".py",
385 "mimetype": "text/x-python",
386 "name": "python",
387 "nbconvert_exporter": "python",
388 "pygments_lexer": "ipython3",
389 "version": "3.8.5"
390 }
391 },
392 "nbformat": 4,
393 "nbformat_minor": 4
394}
diff --git a/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-notebook-checkpoint.ipynb b/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-notebook-checkpoint.ipynb
new file mode 100644
index 0000000..16a6d40
--- /dev/null
+++ b/src/Lecture7/slides/.ipynb_checkpoints/X1-ComputationalComplexity-notebook-checkpoint.ipynb
@@ -0,0 +1,412 @@
1{
2 "cells": [
3 {
4 "cell_type": "markdown",
5 "metadata": {},
6 "source": [
7 "# Nested loops\n",
8 "\n",
9 "The following two functions compute sum and product of matrices, respectively.\n",
10 "\n",
11 "By counting the nested loops it is easy to see that `add()` is $O(n^2)$ while `prod()` is $O(n^3)$."
12 ]
13 },
14 {
15 "cell_type": "code",
16 "execution_count": 2,
17 "metadata": {},
18 "outputs": [
19 {
20 "name": "stdout",
21 "output_type": "stream",
22 "text": [
23 "Time for add: 0.005766554000000035\n",
24 "Time for prod: 1.3871021639999999\n"
25 ]
26 }
27 ],
28 "source": [
29 "from random import randint\n",
30 "import time\n",
31 "\n",
32 "def add(A, B):\n",
33 " S = [[0] * len(A) for i in range(len(A))]\n",
34 " for i in range(len(A)):\n",
35 " for j in range(len(A)):\n",
36 " S[i][j] = A[i][j] + B[i][j]\n",
37 " return S\n",
38 "\n",
39 "def prod(A, B):\n",
40 " S = [[0] * len(A) for i in range(len(A))]\n",
41 " for i in range(len(A)):\n",
42 " for j in range(len(A)):\n",
43 " for k in range(len(A)):\n",
44 " S[i][j] = S[i][j] + A[i][k] * B[k][j]\n",
45 " return S\n",
46 "\n",
47 "N = 200\n",
48 "A = [ [randint(0,100) for i in range(N)] for j in range(N) ]\n",
49 "B = [ [randint(0,100) for i in range(N)] for j in range(N) ]\n",
50 "\n",
51 "t0 = time.process_time()\n",
52 "add(A,B)\n",
53 "t1 = time.process_time()\n",
54 "prod(A,B)\n",
55 "t2 = time.process_time()\n",
56 "\n",
57 "print(\"Time for add: \", t1-t0)\n",
58 "print(\"Time for prod:\", t2-t1)"
59 ]
60 },
61 {
62 "cell_type": "markdown",
63 "metadata": {},
64 "source": [
65 "# Sorting a list, slow version\n",
66 "\n",
67 "The following code implements a slow version of the so-called *insertion sort* alogithm\n",
68 "\n",
69 "Complexity: $O(n^2)$."
70 ]
71 },
72 {
73 "cell_type": "code",
74 "execution_count": 1,
75 "metadata": {},
76 "outputs": [
77 {
78 "name": "stdout",
79 "output_type": "stream",
80 "text": [
81 "Running time: 1.1191449070000001\n"
82 ]
83 }
84 ],
85 "source": [
86 "from random import randint\n",
87 "import time\n",
88 "\n",
89 "def correct_position(e, S):\n",
90 " for i in range(len(S)):\n",
91 " if S[i] > e:\n",
92 " return i\n",
93 " return len(S)\n",
94 "\n",
95 "def sort_list(L):\n",
96 " S = []\n",
97 " for e in L:\n",
98 " cp = correct_position(e, S)\n",
99 " S.insert(cp, e)\n",
100 " return S\n",
101 "\n",
102 "N = 10000\n",
103 "L = [randint(0,10**9) for i in range(N)]\n",
104 "\n",
105 "t0 = time.process_time()\n",
106 "sort_list(L)\n",
107 "t1 = time.process_time()\n",
108 "\n",
109 "print(\"Running time:\", t1-t0)"
110 ]
111 },
112 {
113 "cell_type": "markdown",
114 "metadata": {},
115 "source": [
116 "# Binary search\n",
117 "\n",
118 "The following code implements a binary search.\n",
119 "\n",
120 "Complexity: $O(\\log_2(n))$"
121 ]
122 },
123 {
124 "cell_type": "code",
125 "execution_count": 3,
126 "metadata": {},
127 "outputs": [
128 {
129 "name": "stdout",
130 "output_type": "stream",
131 "text": [
132 "The correct position of e = 216197744 in L is:\n",
133 "... 216196218 216197540 e 216198673 216198962 ...\n",
134 "\n",
135 "Time for sorting: 0.26413054400000036\n",
136 "Time for searching: 9.616099999965044e-05\n"
137 ]
138 }
139 ],
140 "source": [
141 "from random import randint\n",
142 "import time\n",
143 "\n",
144 "def binary_search(e, S, start, end):\n",
145 " if start == end:\n",
146 " return start\n",
147 " midpoint = (start+end) // 2\n",
148 " if e < S[midpoint]:\n",
149 " return binary_search(e, S, start, midpoint)\n",
150 " else:\n",
151 " return binary_search(e, S, midpoint+1, end)\n",
152 " \n",
153 "N = 1000000\n",
154 "L = [randint(0,10**9) for i in range(N)]\n",
155 "e = randint(0,10**9)\n",
156 "\n",
157 "t0 = time.process_time()\n",
158 "L.sort() # Using Python's sort()\n",
159 "t1 = time.process_time()\n",
160 "i = binary_search(e, L, 0, len(L))\n",
161 "t2 = time.process_time()\n",
162 "print(\"The correct position of e =\", e, \"in L is:\")\n",
163 "print(\"...\", L[i-2], L[i-1], \"e\", L[i], L[i+1], \"...\")\n",
164 "print(\"\")\n",
165 "print(\"Time for sorting: \", t1-t0)\n",
166 "print(\"Time for searching:\", t2-t1)\n"
167 ]
168 },
169 {
170 "cell_type": "markdown",
171 "metadata": {},
172 "source": [
173 "# Sorting a list, fast version (with binary_search)\n",
174 "\n",
175 "The following code uses the function `binary_search()` above instead of `correct_position()` in our insertion sort algorithm.\n",
176 "\n",
177 "Complexity: $O(n\\log_2(n))$"
178 ]
179 },
180 {
181 "cell_type": "code",
182 "execution_count": 69,
183 "metadata": {},
184 "outputs": [
185 {
186 "name": "stdout",
187 "output_type": "stream",
188 "text": [
189 "Running time: 0.03710268399998995\n"
190 ]
191 }
192 ],
193 "source": [
194 "from random import randint\n",
195 "import time\n",
196 "\n",
197 "def binary_search(e, S, start, end):\n",
198 " if start == end:\n",
199 " return start\n",
200 " midpoint = (start+end) // 2\n",
201 " if e < S[midpoint]:\n",
202 " return binary_search(e, S, start, midpoint)\n",
203 " else:\n",
204 " return binary_search(e, S, midpoint+1, end)\n",
205 " \n",
206 "def sort_list(L):\n",
207 " S = []\n",
208 " for e in L:\n",
209 " cp = binary_search(e, S, 0, len(S)) # Changed here\n",
210 " S.insert(cp, e)\n",
211 " return S\n",
212 " \n",
213 "N = 10000\n",
214 "L = [randint(0,10**9) for i in range(N)]\n",
215 "\n",
216 "t0 = time.process_time()\n",
217 "sort_list(L)\n",
218 "t1 = time.process_time()\n",
219 "\n",
220 "print(\"Running time:\", t1-t0)"
221 ]
222 },
223 {
224 "cell_type": "markdown",
225 "metadata": {},
226 "source": [
227 "# Fast exponentiation\n",
228 "\n",
229 "The following cell contains two functions for computing $a^n$ ($n$ non-negative integer): a slow one that runs in $O(n)$ and a fast one that runs in $O(\\log_2(n))$. We compare these two also with Python's built-in operator `**`.\n",
230 "\n",
231 "Complexity: $O(n)$ for the slow algorithm, $O(\\log_2(n))$ for the other two."
232 ]
233 },
234 {
235 "cell_type": "code",
236 "execution_count": 30,
237 "metadata": {},
238 "outputs": [
239 {
240 "name": "stdout",
241 "output_type": "stream",
242 "text": [
243 "2.71828179834636\n",
244 "2.7182817863957984\n",
245 "2.7182817983473577\n",
246 "Time for slow_power(): 3.234879998000004\n",
247 "Time for fast_power(): 9.059099999575437e-05\n",
248 "Time for Python's **: 0.00010159500000384014\n"
249 ]
250 }
251 ],
252 "source": [
253 "import time\n",
254 "\n",
255 "def slow_power(a, n):\n",
256 " r = 1\n",
257 " for i in range(n):\n",
258 " r = r * a\n",
259 " return r\n",
260 "\n",
261 "def fast_power(a, n):\n",
262 " if n == 0:\n",
263 " return 1\n",
264 " if n%2 == 0:\n",
265 " return fast_power(a*a, n//2)\n",
266 " else:\n",
267 " return a * fast_power(a, n-1)\n",
268 "\n",
269 "a = 1.00000001\n",
270 "n = 100000000\n",
271 "\n",
272 "t0 = time.process_time()\n",
273 "print(slow_power(a, n))\n",
274 "t1 = time.process_time()\n",
275 "print(fast_power(a, n))\n",
276 "t2 = time.process_time()\n",
277 "print(a**n)\n",
278 "t3 = time.process_time()\n",
279 "\n",
280 "print(\"Time for slow_power():\", t1-t0)\n",
281 "print(\"Time for fast_power():\", t2-t1)\n",
282 "print(\"Time for Python's **: \", t3-t2)"
283 ]
284 },
285 {
286 "cell_type": "markdown",
287 "metadata": {},
288 "source": [
289 "# Fast gcd\n",
290 "\n",
291 "Complexity: $O(\\log_2(n))$"
292 ]
293 },
294 {
295 "cell_type": "code",
296 "execution_count": 31,
297 "metadata": {},
298 "outputs": [
299 {
300 "name": "stdout",
301 "output_type": "stream",
302 "text": [
303 "126\n",
304 "Running time: 0.00017707599999994272\n"
305 ]
306 }
307 ],
308 "source": [
309 "import time\n",
310 "\n",
311 "def gcd(a, b):\n",
312 " if b == 0:\n",
313 " return a\n",
314 " else:\n",
315 " return gcd(b, a%b)\n",
316 "\n",
317 "t0 = time.process_time()\n",
318 "print(gcd(155275387236018, 572335397352432))\n",
319 "t1 = time.process_time()\n",
320 "\n",
321 "print(\"Running time:\", t1-t0)"
322 ]
323 },
324 {
325 "cell_type": "markdown",
326 "metadata": {},
327 "source": [
328 "# Fibonacci numbers\n",
329 "\n",
330 "In the following cell there are two functions that compute the $n$-th Fibonacci number. They are almost the same, but the second one memorizes the results in a list to avoid computing them multiple times, and it is much much faster.\n",
331 "\n",
332 "Complexity: $O\\left(\\left(\\frac{1+\\sqrt 5}{2}\\right)^n\\right)\\sim O(1.6^n)$ for the slow version, $O(n)$ for the fast version."
333 ]
334 },
335 {
336 "cell_type": "code",
337 "execution_count": 37,
338 "metadata": {},
339 "outputs": [
340 {
341 "name": "stdout",
342 "output_type": "stream",
343 "text": [
344 "9227465\n",
345 "9227465\n",
346 "Time for F_slow: 2.3301570169999906\n",
347 "Time for F_fast: 8.848800000293977e-05\n"
348 ]
349 }
350 ],
351 "source": [
352 "import time\n",
353 "\n",
354 "F_memorized = [-1] * (10**6)\n",
355 "\n",
356 "def F_slow(n):\n",
357 " if n <= 1:\n",
358 " return n\n",
359 " else:\n",
360 " return F_slow(n-1) + F_slow(n-2)\n",
361 " \n",
362 "def F_fast(n):\n",
363 " if F_memorized[n] == -1:\n",
364 " if n <= 1:\n",
365 " F_memorized[n] = n\n",
366 " else:\n",
367 " F_memorized[n] = F_fast(n-1) + F_fast(n-2)\n",
368 " \n",
369 " return F_memorized[n]\n",
370 "\n",
371 "n = 35\n",
372 "\n",
373 "t0 = time.process_time()\n",
374 "print(F_slow(n))\n",
375 "t1 = time.process_time()\n",
376 "print(F_fast(n))\n",
377 "t2 = time.process_time()\n",
378 "\n",
379 "print(\"Time for F_slow:\", t1-t0)\n",
380 "print(\"Time for F_fast:\", t2-t1)"
381 ]
382 },
383 {
384 "cell_type": "code",
385 "execution_count": null,
386 "metadata": {},
387 "outputs": [],
388 "source": []
389 }
390 ],
391 "metadata": {
392 "kernelspec": {
393 "display_name": "Python 3",
394 "language": "python",
395 "name": "python3"
396 },
397 "language_info": {
398 "codemirror_mode": {
399 "name": "ipython",
400 "version": 3
401 },
402 "file_extension": ".py",
403 "mimetype": "text/x-python",
404 "name": "python",
405 "nbconvert_exporter": "python",
406 "pygments_lexer": "ipython3",
407 "version": "3.8.5"
408 }
409 },
410 "nbformat": 4,
411 "nbformat_minor": 4
412}
diff --git a/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-checkpoint.ipynb b/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-checkpoint.ipynb
new file mode 100644
index 0000000..a35bb8d
--- /dev/null
+++ b/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-checkpoint.ipynb
@@ -0,0 +1,165 @@
1{
2 "cells": [
3 {
4 "cell_type": "markdown",
5 "metadata": {},
6 "source": [
7 "# Diffie-Hellman key exchange\n",
8 "\n",
9 "The following is a simple implementation of the classic [Diffie-Hellman key exchange](https://en.wikipedia.org/wiki/Diffie%E2%80%93Hellman_key_exchange) cryptographic protocol."
10 ]
11 },
12 {
13 "cell_type": "code",
14 "execution_count": 9,
15 "metadata": {},
16 "outputs": [
17 {
18 "name": "stdout",
19 "output_type": "stream",
20 "text": [
21 "Public key: p = 20747 and g = 13428 \n",
22 "\n",
23 "[[ Alice's secret key: a = 12403 ]]\n",
24 "[[ Bob's secret key: b = 17642 ]] \n",
25 "\n",
26 "Alice sends h1 = 14710 to Bob\n",
27 "Bob sends h2 = 10680 to Alice \n",
28 "\n",
29 "Alice computed 10455 using h2 and her secret a\n",
30 "Bob computed 10455 using h1 and his secret b\n"
31 ]
32 }
33 ],
34 "source": [
35 "# Public information:\n",
36 "p = Primes()[10^3 + randint(1,10000)] # random prime\n",
37 "g = randint(2, p-1) # random integer\n",
38 "\n",
39 "print(\"Public key: p =\", p, \"and g =\", g, \"\\n\")\n",
40 "\n",
41 "a = randint(2, p-1) # Only Alice knows this\n",
42 "b = randint(2, p-1) # Only Bob knows this\n",
43 "\n",
44 "print(\"[[ Alice's secret key: a =\", a, \"]]\")\n",
45 "print(\"[[ Bob's secret key: b =\", b, \"]]\", \"\\n\")\n",
46 "\n",
47 "h1 = (g^a) % p # Alice sends this to Bob\n",
48 "h2 = (g^b) % p # Bob sends this to Alice\n",
49 "\n",
50 "print(\"Alice sends h1 =\", h1, \"to Bob\")\n",
51 "print(\"Bob sends h2 =\", h2, \"to Alice\", \"\\n\")\n",
52 "\n",
53 "secret_a = (h2^a) % p # Alice can compute this because she knows a\n",
54 "secret_b = (h1^b) % p # Bob can compute this because he knows b\n",
55 "\n",
56 "print(\"Alice computed\", secret_a, \"using h2 and her secret a\")\n",
57 "print(\"Bob computed\", secret_b, \"using h1 and his secret b\")"
58 ]
59 },
60 {
61 "cell_type": "markdown",
62 "metadata": {},
63 "source": [
64 "## General Diffie-Hellman\n",
65 "\n",
66 "The following code is an implementation of a generic Diffie-Hellman key exchange protocol that uses a group $G$ instead of $(\\mathbb Z/p \\mathbb Z)^\\times$."
67 ]
68 },
69 {
70 "cell_type": "code",
71 "execution_count": 24,
72 "metadata": {},
73 "outputs": [
74 {
75 "name": "stdout",
76 "output_type": "stream",
77 "text": [
78 "Public key:\n",
79 "G = Additive abelian group isomorphic to Z/171 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x + 156 over Finite Field of size 157 \n",
80 "g = (155 : 60 : 1) \n",
81 "\n",
82 "[[ Alice's secret key: a = 141 ]]\n",
83 "[[ Bob's secret key: b = 158 ]] \n",
84 "\n",
85 "Alice sends h1 = (29 : 125 : 1) to Bob\n",
86 "Bob sends h2 = (60 : 59 : 1) to Alice \n",
87 "\n",
88 "Alice computed (109 : 94 : 1) using h2 and her secret a\n",
89 "Bob computed (109 : 94 : 1) using h1 and his secret b\n"
90 ]
91 }
92 ],
93 "source": [
94 "def genericDH(G):\n",
95 " if G.cardinality() == 1:\n",
96 " print(\"Group is trivial, can't do anything\")\n",
97 " return\n",
98 " g = G.random_element()\n",
99 " while g == G.identity(): # Make sure g is not trivial\n",
100 " g = G.random_element()\n",
101 " \n",
102 " print(\"Public key:\\nG =\", G, \"\\ng =\", g, \"\\n\")\n",
103 " \n",
104 " a = randint(2, G.exponent()-1) # Only Alice knows this\n",
105 " b = randint(2, G.exponent()-1) # Only Bob knows this\n",
106 "\n",
107 " print(\"[[ Alice's secret key: a =\", a, \"]]\")\n",
108 " print(\"[[ Bob's secret key: b =\", b, \"]]\", \"\\n\")\n",
109 " \n",
110 " # \"Ternary operator\", I did not explain this\n",
111 " # https://docs.python.org/3/reference/expressions.html#conditional-expressions\n",
112 " h1 = g^a if G.is_multiplicative() else a*g # Alice sends this to Bob\n",
113 " h2 = g^b if G.is_multiplicative() else b*g # Bob sends this to Alice\n",
114 "\n",
115 " print(\"Alice sends h1 =\", h1, \"to Bob\")\n",
116 " print(\"Bob sends h2 =\", h2, \"to Alice\", \"\\n\")\n",
117 " \n",
118 " secret_a = h2^a if G.is_multiplicative() else a*h2 # Alice can compute this because she knows a\n",
119 " secret_b = h1^b if G.is_multiplicative() else b*h1 # Bob can compute this because he knows b\n",
120 "\n",
121 " print(\"Alice computed\", secret_a, \"using h2 and her secret a\")\n",
122 " print(\"Bob computed\", secret_b, \"using h1 and his secret b\")\n",
123 " \n",
124 "E = EllipticCurve(GF(157), [1,-1])\n",
125 "G = E.abelian_group()\n",
126 "genericDH(G)"
127 ]
128 },
129 {
130 "cell_type": "markdown",
131 "metadata": {},
132 "source": [
133 "# Numerical methods for PDEs"
134 ]
135 },
136 {
137 "cell_type": "code",
138 "execution_count": null,
139 "metadata": {},
140 "outputs": [],
141 "source": []
142 }
143 ],
144 "metadata": {
145 "kernelspec": {
146 "display_name": "SageMath 9.0",
147 "language": "sage",
148 "name": "sagemath"
149 },
150 "language_info": {
151 "codemirror_mode": {
152 "name": "ipython",
153 "version": 3
154 },
155 "file_extension": ".py",
156 "mimetype": "text/x-python",
157 "name": "python",
158 "nbconvert_exporter": "python",
159 "pygments_lexer": "ipython3",
160 "version": "3.8.5"
161 }
162 },
163 "nbformat": 4,
164 "nbformat_minor": 4
165}
diff --git a/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-notebook-checkpoint.ipynb b/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-notebook-checkpoint.ipynb
new file mode 100644
index 0000000..ccd6a21
--- /dev/null
+++ b/src/Lecture7/slides/.ipynb_checkpoints/X2-StudentsRequests-notebook-checkpoint.ipynb
@@ -0,0 +1,298 @@
1{
2 "cells": [
3 {
4 "cell_type": "markdown",
5 "metadata": {},
6 "source": [
7 "# Diffie-Hellman key exchange\n",
8 "\n",
9 "The following is a simple implementation of the classic [Diffie-Hellman key exchange](https://en.wikipedia.org/wiki/Diffie%E2%80%93Hellman_key_exchange) cryptographic protocol."
10 ]
11 },
12 {
13 "cell_type": "code",
14 "execution_count": 1,
15 "metadata": {},
16 "outputs": [
17 {
18 "name": "stdout",
19 "output_type": "stream",
20 "text": [
21 "Public key: p = 75521 and g = 58258 \n",
22 "\n",
23 "[[ Alice's secret key: a = 22794 ]]\n",
24 "[[ Bob's secret key: b = 69773 ]] \n",
25 "\n",
26 "Alice sends h1 = 31067 to Bob\n",
27 "Bob sends h2 = 54398 to Alice \n",
28 "\n",
29 "Alice computed 30031 using h2 and her secret a\n",
30 "Bob computed 30031 using h1 and his secret b\n"
31 ]
32 }
33 ],
34 "source": [
35 "# Public information:\n",
36 "p = Primes()[10^3 + randint(1,10000)] # random prime\n",
37 "g = randint(2, p-1) # random integer\n",
38 "\n",
39 "print(\"Public key: p =\", p, \"and g =\", g, \"\\n\")\n",
40 "\n",
41 "a = randint(2, p-1) # Only Alice knows this\n",
42 "b = randint(2, p-1) # Only Bob knows this\n",
43 "\n",
44 "print(\"[[ Alice's secret key: a =\", a, \"]]\")\n",
45 "print(\"[[ Bob's secret key: b =\", b, \"]]\", \"\\n\")\n",
46 "\n",
47 "h1 = (g^a) % p # Alice sends this to Bob\n",
48 "h2 = (g^b) % p # Bob sends this to Alice\n",
49 "\n",
50 "print(\"Alice sends h1 =\", h1, \"to Bob\")\n",
51 "print(\"Bob sends h2 =\", h2, \"to Alice\", \"\\n\")\n",
52 "\n",
53 "secret_a = (h2^a) % p # Alice can compute this because she knows a\n",
54 "secret_b = (h1^b) % p # Bob can compute this because he knows b\n",
55 "\n",
56 "print(\"Alice computed\", secret_a, \"using h2 and her secret a\")\n",
57 "print(\"Bob computed\", secret_b, \"using h1 and his secret b\")"
58 ]
59 },
60 {
61 "cell_type": "markdown",
62 "metadata": {},
63 "source": [
64 "## General Diffie-Hellman\n",
65 "\n",
66 "The following code is an implementation of a generic Diffie-Hellman key exchange protocol that uses a group $G$ instead of $(\\mathbb Z/p \\mathbb Z)^\\times$."
67 ]
68 },
69 {
70 "cell_type": "code",
71 "execution_count": 2,
72 "metadata": {},
73 "outputs": [
74 {
75 "name": "stdout",
76 "output_type": "stream",
77 "text": [
78 "Public key:\n",
79 "G = Additive abelian group isomorphic to Z/171 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x + 156 over Finite Field of size 157 \n",
80 "g = (53 : 90 : 1) \n",
81 "\n",
82 "[[ Alice's secret key: a = 145 ]]\n",
83 "[[ Bob's secret key: b = 65 ]] \n",
84 "\n",
85 "Alice sends h1 = (150 : 80 : 1) to Bob\n",
86 "Bob sends h2 = (4 : 58 : 1) to Alice \n",
87 "\n",
88 "Alice computed (28 : 28 : 1) using h2 and her secret a\n",
89 "Bob computed (28 : 28 : 1) using h1 and his secret b\n"
90 ]
91 }
92 ],
93 "source": [
94 "def genericDH(G):\n",
95 " if G.cardinality() == 1:\n",
96 " print(\"Group is trivial, can't do anything\")\n",
97 " return\n",
98 " g = G.random_element()\n",
99 " while g == G.identity(): # Make sure g is not trivial\n",
100 " g = G.random_element()\n",
101 " \n",
102 " print(\"Public key:\\nG =\", G, \"\\ng =\", g, \"\\n\")\n",
103 " \n",
104 " a = randint(2, G.exponent()-1) # Only Alice knows this\n",
105 " b = randint(2, G.exponent()-1) # Only Bob knows this\n",
106 "\n",
107 " print(\"[[ Alice's secret key: a =\", a, \"]]\")\n",
108 " print(\"[[ Bob's secret key: b =\", b, \"]]\", \"\\n\")\n",
109 " \n",
110 " # \"Ternary operator\", I did not explain this\n",
111 " # https://docs.python.org/3/reference/expressions.html#conditional-expressions\n",
112 " h1 = g^a if G.is_multiplicative() else a*g # Alice sends this to Bob\n",
113 " h2 = g^b if G.is_multiplicative() else b*g # Bob sends this to Alice\n",
114 "\n",
115 " print(\"Alice sends h1 =\", h1, \"to Bob\")\n",
116 " print(\"Bob sends h2 =\", h2, \"to Alice\", \"\\n\")\n",
117 " \n",
118 " secret_a = h2^a if G.is_multiplicative() else a*h2 # Alice can compute this because she knows a\n",
119 " secret_b = h1^b if G.is_multiplicative() else b*h1 # Bob can compute this because he knows b\n",
120 "\n",
121 " print(\"Alice computed\", secret_a, \"using h2 and her secret a\")\n",
122 " print(\"Bob computed\", secret_b, \"using h1 and his secret b\")\n",
123 " \n",
124 "E = EllipticCurve(GF(157), [1,-1])\n",
125 "G = E.abelian_group()\n",
126 "genericDH(G)"
127 ]
128 },
129 {
130 "cell_type": "markdown",
131 "metadata": {},
132 "source": [
133 "# Numerical methods for differential equations\n",
134 "\n",
135 "## Euler's method (ODE)\n",
136 "\n",
137 "In sage you can use [`ode_solver()`](https://doc.sagemath.org/html/en/reference/calculus/sage/calculus/ode.html) to solve any ordinary differential equation by hand, but Euler's method is very simple to implement by hand:"
138 ]
139 },
140 {
141 "cell_type": "code",
142 "execution_count": 4,
143 "metadata": {},
144 "outputs": [
145 {
146 "data": {
147 "image/png": "iVBORw0KGgoAAAANSUhEUgAAAk4AAAGGCAYAAACNCg6xAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3de5yOdf7H8fcYzFCMs3EYzCLHQowIlYTVQee2ItRurUUl20ltRduGTmtpapeKWp1+5dhBpeQUyjGHEFKImRznhBlzz/X747NjHIa5cc993fd9vZ6Px/24m9M1n7lNeff9fq7PN8pxHEcAAAAoUgm3CwAAAAgXBCcAAAA/EZwAAAD8RHACAADwE8EJAADATwQnAAAAPxGcAAAA/ERwAgAA8BPBCfCT4zhKT08XM2MBwLsIToCfMjIyFBcXp4yMDLdLAQC4hOAEAADgJ4ITAACAnwhOAAAAfiI4AQAA+IngBAAAIs6XX0qzZwf+ugQnoAjJyclq2rSpkpKS3C4FAOCHnBypf3/phRcCf+0oh6E0gF/S09MVFxentLQ0lS9f3u1yAAAnMXasNHiw9P33UvPmgb02K04AACBipKVJTz8t9esX+NAkEZwAAEAEGTVKysqy8FQcCE4AACAibN8u/fOf0pAhUq1axfM9CE4AACAiPPmkVK6c9PDDxfc9ShbfpQEAAIJj1Spp4kRrDC/O+3e4qw7wE3fVAUBochypa1dp2zZpzRqpVKni+16sOAEAgLD26afSV19J06cXb2iS6HECisQATAAIXYcPSw8+KF1+uXTNNcX//diqA/zEVh0AhJ6XX5buu09avlxq2bL4vx8rTgAAICzt2ycNGybdeWdwQpNEcAIAAGHqH/+QDh2SnnkmeN+T4AQAAMLO5s3SmDHSo49KNWoE7/sSnOAJ9erVU1RU1AmPgQMHul0aAOAMPPKIVL26TQkPJsYRwBOWLFkin8935O01a9aoa9euuvnmm12sCgBwJubNkyZPliZNksqWDe735q46eNLgwYP18ccfa+PGjYqKivLra7irDgDcl5cntW0rlSghLV5sz8HEihM8JycnR5MmTdKQIUNOGZqys7OVnZ195O309PRglAcAOIW335aWLZPmzw9+aJLocYIHTZs2Tfv371e/fv1O+XkjRoxQXFzckUdCQkJwCgQAFCojw3qbbr5Z6tjRnRrYqoPndO/eXaVLl9ZHH310ys8rbMUpISGBrToAcMnQodLo0dL69VLduu7UwFYdPOWXX37Rl19+qSlTphT5uTExMYqJiQlCVQCAomzcKL30kvTYY+6FJomtOnjMhAkTVK1aNV111VVulwIAOA1Dhti8pocfdrcOVpzgGXl5eZowYYL69u2rkiX51QeAcDFzpvTxx9KHH0plyrhbCytO8Iwvv/xSW7du1V133eV2KQAAP+XkSIMHS5dfLt1wg9vVsOIED+nWrZu4FwIAwsuYMXa8yuTJkp9j94oVK04AACAk7dwpDR8uDRggNW/udjWG4AQAAELS0KFSTIyFp1BBcAKKkJycrKZNmyopKcntUgDAM779VnrzTenZZ6WKFd2upgADMAE/cVYdAARHXp7Urp2UmystWSJFR7tdUQGawwEAQEh5800LTPPnh1ZoktiqAwAAISQtTXr0Uen22907j+5UCE4AACBkPPGElJUljRrldiWFY6sOAACEhBUrpORk6bnnpNq13a6mcDSHA36iORwAik9ennTxxVJmpgWoUqXcrqhwrDgBAADXvfGGjSCYNy90Q5NEjxMAAHDZ7t3SI49IffpInTq5Xc2pEZyAIjAAEwCK19Chks9nvU2hjh4nwE/0OAFA4C1aZL1Nycl2Jl2oIzgBfiI4AUBg5eZKSUk25PLbb0Nv2GVhaA4HAACueOUV6fvvpcWLwyM0SfQ4AQAAF+zYYcMu77lHatvW7Wr8R3ACAABBd999Umys9OyzbldyetiqAwAAQfXRR9LkydK770qVKrldzemhORzwE83hAHD2MjKkZs3s8emnUlSU2xWdHrbqAABA0DzxhLRnjzWGh1tokghOQJEYgAkAgbFkiTR2rDR8uJSY6HY1Z4atOsBPbNUBwJnLn9kkWYAqGaZd1mFaNgAACCejR0urVtnMpnANTRJbdQAAoJht2SI99ZR0770Fq07hiuAEAACKjePYGXSVK0t//7vb1Zy9MF4sAwAAoe7996XPPpNmzJDKlXO7mrPHihMAACgW+/ZJ998v3XSTdM01blcTGAQnAABQLB5+WDp0SPrXv9yuJHDYqgMAAAH31VfSa69Jr74q1azpdjWBw4oTUAQGYALA6cnMlO6+W7rsMumee9yuJrAYgAn4iQGYAOCf+++Xxo+XVq+W6td3u5rAYqsOAAAEzIIFdqzKiy9GXmiSWHEC/MaKEwCc2sGDUsuWUqVKFqCio92uKPBYcQIAAAExfLj088/StGmRGZokmsPhIb/++qt69+6typUrq2zZsmrZsqWWLVvmdlkAEBGWLJGef14aNkxq0sTtaooPK07whH379qlDhw7q3LmzZs6cqWrVqmnz5s2qUKGC26UBQNjLyZHuuktq0UJ68EG3qyleBCd4wqhRo5SQkKAJEyYceV+9evXcKwgAIsizz0rr19uqU6lSbldTvNiqgyfMmDFDbdq00c0336xq1aqpVatWGj9+vNtlAUDYW7VK+sc/pKFDrTE80nFXHTwhNjZWkjRkyBDdfPPN+u677zR48GD95z//UZ8+fQr9muzsbGVnZx95Oz09XQkJCdxVBwD/k5MjtW0r+XzS0qVSTIzbFRU/ghM8oXTp0mrTpo0WLlx45H333XeflixZokWLFhX6NcOGDdPw4cNPeD/BCQDME09II0dK330ntWrldjXBwVYdPKFGjRpq2rTpMe9r0qSJtm7detKvGTp0qNLS0o48tm3bVtxlAkDY+O47acQI6cknvROaJJrD4REdOnTQhg0bjnnfjz/+qLp16570a2JiYhTjhXVnADhNBw9KfftaYHr0UberCS6CEzzhgQce0MUXX6xnn31Wt9xyi7777juNGzdO48aNc7s0AAg7f/ubtGWLtHx55N9Fdzx6nOAZH3/8sYYOHaqNGzcqMTFRQ4YM0d133+3313PkCgBI8+dLl14qPfdc5M9sKgzBCfATwQmA12Vm2pDLGjWkuXMj91iVU2GrDgAA+OXhh6WUFOnzz70ZmiSCEwAA8MOsWdKrr0rJyVKDBm5X4x626gA/sVUHwKv27ZMuuEBq3NhWm0p4eJiRh390wD/Jyclq2rSpkpKS3C4FAILOcaT+/a2/6Y03vB2aJFacAL+x4gTAi956y2Y2vf++dMstblfjPo/nRgAAcDJbtkiDBkl9+hCa8hGcAADACXJzpTvukCpXlsaOdbua0MFddQAA4AQjR0qLFknz5kl0JxRgxQkAABzju++kYcOkxx6TOnRwu5rQQnM44CeawwF4QWamHd5bsaL0zTfeO4uuKGzVAQCAIx54QNqxQ/rkE0JTYQhOAABAkjR1qvTaa9K4cdJ557ldTWiixwkoAgMwAXjBzp3S3XdL114r/elPblcTuuhxAvxEjxOASJWXJ3XrJq1dK61aJVWt6nZFoYutOgAAPG7UKGn2bOmLLwhNRWGrDgAAD1u4UHriCWnoUOmKK9yuJvSxVQf4ia06AJFm3z6pZUupdm1pzhzuovMHK04AAHiQ41gTeHq69M47hCZ/0eMEAIAH/fvf0pQp0uTJUt26blcTPlhxAgDAY1atskGXAwZIN9zgdjXhhR4nwE/0OAGIBFlZUps2UunS0rffSrGxblcUXtiqA4qQnJys5ORk+Xw+t0sBgLN2333S1q3S0qWEpjPBihPgJ1acAIS7t9+WeveW3nhDuvNOt6sJT/Q4AQDgAT/8IN1zjwWnfv3criZ8seIE+IkVJwDhKjNTattWKlHC+prOOcftisIXPU4AAEQwx5H697e+piVLCE1ni+AEAEAEGzfOepvefVdq0sTtasIfPU4AAESoZcvsLroBA6Rbb3W7mshAjxPgJ3qcAISTffuk1q2lypWlBQukmBi3K4oMbNUBABBhHMfGDezbJ331FaEpkNiqA4qQnJyspk2bKikpye1SAMAvL74oTZ8uvfWWlJjodjWRha06wE9s1QEIB3PnSl26SA8+KI0c6XY1kYfgBPiJ4AQg1G3fbn1NzZtLn38ulaQhJ+DYqgMAIAJkZ0s33mjnz733HqGpuPCyAgAQAe69V/r+e7uDrmpVt6uJXKw4wROGDRumqKioYx7x8fFulwUAATF+vD1efVVq08btaiIbK07wjGbNmunLL7888nZ0dLSL1QBAYHz7rTRokPSXv9gIAhQvghM8o2TJkqwyAYgoqanW19S6tTR6tNvVeANbdfCMjRs3qmbNmkpMTNStt96qn3766ZSfn52drfT09GMeABAqDh+WbrlF8vmkDz+USpd2uyJvIDjBEy666CK99dZb+vzzzzV+/HilpKTo4osv1p49e076NSNGjFBcXNyRR0JCQhArBoBTe+ghaeFC6YMPpJo13a7GO5jjBE/KyspS/fr19fDDD2vIkCGFfk52drays7OPvJ2enq6EhATmOAFw3ZtvSv36SWPHWn8TgoceJ3jSOeeco/PPP18bN2486efExMQohgOeAISYRYuke+6R/vhHaeBAt6vxHrbq4EnZ2dlat26datSo4XYpAOC37dul66+X2raVXnlFiopyuyLvITjBEx588EHNnTtXW7Zs0bfffqubbrpJ6enp6tu3r9ulAYBfDhyQrr1WiomRJk+mGdwtbNXBE7Zv367bbrtNu3fvVtWqVdWuXTstXrxYdevWdbs0ACiS40h33SWtXy99841UrZrbFXkXwQme8N5777ldAgCcsWefld5/38YOtGzpdjXexlYdAAAhbPp06W9/k556yoZdwl0EJwAAQtTq1VKvXhaYnnzS7WogEZyAIiUnJ6tp06ZKSkpyuxQAHpKaKvXsKTVoYHObSvA3dkhgACbgp/T0dMXFxTEAE0CxO3hQ6txZ+uUXO8S3Th23K0I+msMBAAgheXlS377SqlXSvHmEplBDcAIAIIQ88YTdPTd5stSmjdvV4HgEJwAAQsTEiTZ64LnnbEI4Qg+tZgAAhIA5c+wMurvvlh580O1qcDIEJwAAXLZhg3TDDdKll0rJyZxBF8oITgAAuGj3bumqq6QaNaQPPpBKlXK7IpwKPU4AALjk0CHrZUpPt7EDFSq4XRGKwooTUAQGYAIoDnl5Up8+0tKl0owZUmKi2xXBHwzABPzEAEwAgfTAA9KYMTZ24Lrr3K4G/mKrDgCAIHvpJWn0aGsEJzSFF7bqAAAIovffl/76V+nRR6UBA9yuBqeL4AQAQJDMmWN9Tb1726BLhB+CEwAAQbBmjW3Ldeokvf46s5rCFcEJAIBitn271KOHVLeuNYOXLu12RThTBCcAAIrRnj1St25SiRLSp59KcXFuV4SzwV11AAAUk8xMmwq+a5e0YIFUq5bbFeFsseIEFIEBmADORE6OdOON0tq10mefSY0auV0RAoEBmICfGIAJwF8+n9SrlzR1qjRzpnT55W5XhEBhqw4AgAByHOm+++zA3g8+IDRFGoITAAABNGyY9Mor0vjx0g03uF0NAo0eJwAAAmTsWOnpp6URI6Q//cntalAcCE4AAATApEm2RffXv0qPPOJ2NSguBCcAAM7SlClSv37SnXdKzz3HVPBIRnACAOAsfPqpdOut0s03W19TCf5mjWj88QIAcIZmz7YG8CuvlN56S4qOdrsiFDeCE1AEBmACKMw330g9e0qXXiq9/75UqpTbFSEYGIAJ+IkBmADyLV0qdekitWxpAy7LlnW7IgQLK04AAJyG1aul7t2lJk2kjz8mNHkNwQkAAD9t2CB17SrVqWPnz5Ur53ZFCDaCEwAAftiyxbbnKleWvvhCqlDB7YrgBoITPGnEiBGKiorS4MGD3S4FQBj4+Wepc2epTBnpyy+lqlXdrghuITjBc5YsWaJx48bpggsucLsUAGHg55+lyy6TSpa08QM1arhdEdxEcIKnZGZmqlevXho/frwqVqzodjkAQtyWLTZuoGRJac4cKSHB7YrgNoITPGXgwIG66qqrdMUVV7hdCoAQt2WLrTSVLm2hqXZttytCKCjpdgFAsLz33ntavny5lixZ4tfnZ2dnKzs7+8jb6enpxVUagBDz00/W01S6tPT114QmFGDFCZ6wbds23X///Zo0aZJiY2P9+poRI0YoLi7uyCOBNXrAE376iZUmnByTw+EJ06ZN0/XXX6/oow6S8vl8ioqKUokSJZSdnX3Mx6TCV5wSEhKYHA5EsM2bbaUpNtZWmmrVcrsihBq26uAJXbp00erVq49535133qnGjRvrkUceOSE0SVJMTIxiYmKCVSIAlxGa4A+CEzyhXLlyat68+THvO+ecc1S5cuUT3g/AezZskK64wuY0zZkj1azpdkUIVfQ4AQA8bdUq6ZJLpPLlCU0RZcYM6YEH7DmA6HEC/JSenq64uDh6nIAI8t130u9/L9WrZ8eoVKnidkUIiGnTpOuvl0qUkPLypOnTpZ49A3JpVpwAAJ40b56dPdekiU0EJzRFiCVLpP797Z/z8qToaFtKDBCCEwDAcz7/3Faa2ra1f+bA3giwZ48Fposusg5/yUKTz2fzJQKE4AQA8JSpU6VrrrHVpk8+kc491+2KcFby8qTXXpMaNZLefVcaPVratMm25+67L6DbdBI9TkCRkpOTlZycLJ/Ppx9//JEeJyCMTZok9esn3Xij/XOpUm5XhLOybJk0YIA1q/XpI40aJcXHF+u3JDgBfqI5HAhv48bZTk6/ftL48baLgzC1d6/0+OPSf/4jNW8uJSdLnToF5VuzVQcAiHgvvCD9+c/SoEG2q0NoClN5edLrr0vnnSe98470z39Ky5cHLTRJBCcAQATLy5MefFB66CFboPjXv+wOdYSh5culiy+W/vQn6corbWrp/fdLJYM7y5tfHwBARDp82LblXnxRGjNGeuYZKSrK7apw2vbutT6mNm2krCxp7lzprbeKvZfpZDhyBQAQcbKypFtukWbNshutbr3V7Ypw2vLypIkTpUcekbKzpZdekgYOdL2jn+AEAIgoe/ZIV18trV5t4wa6dnW7Ipy25cstJC1eLPXqJT3/vFSjhttVSWKrDgAQQbZulTp2tDE+c+YQmsLOvn3WwZ+UJGVk2B/ipEkhE5okVpwAABHihx+k7t2tV/ibb+zGK4SJvDzpzTdtW+7QIbsNctAg17flCsOKE1CE5ORkNW3aVElJSW6XAuAkFi2ylaZKlQhNYWfFCvvDu+suWyJcv1564IGQDE0SAzABvzEAEwhN06ZJt99uuzvTp3PuXNjYv1964gnplVekxo1tiGUAz5QrLqw4AQDC1ujR0g03WDM4h/WGify75c47z56ff15auTIsQpNEcAIAhCGfz85vfeAB6eGHpffek2Jj3a4KRVq50qZ833mndMUVNsRyyJCQ3ZYrDM3hAICwkpUl3Xab9OmndlTZPfe4XRGKtH+/9OSTth3XqJE0e7bUubPbVZ0RghMAIGzs3Cldc40tVHz8sfT737tdEU7JcaT//tfOvDlwQBo1yo5JCaMVpuMRnAAAYWHtWjuizOeTFiyQWrRwuyKc0vff20iBBQtsdPsLL0i1arld1VmjxwkAEPK+/NLOd61Y0YZJE5pCWFqarSpdeKGNcf/qKzv3JgJCk0RwAgCEuAkTpB49LDjNny/Vru12RShU/rZco0bS669LI0daM/jll7tdWUARnIAiMAATcIfPZ4Ok77pL+uMfpY8+ksqVc7sqFGrVKumSS6Q+fWyswPr11tdUurTblQUcAzABPzEAEwie9HQbajlzpvTii7bzExXldlU4QVqa9NRT0ssvSw0b2nOXLm5XVaxoDgcAhJSffrI753791UYOdO/udkU4geNIb78tPfiglJkpjRhh6TYCV5iOx1YdACBkzJkjtW0r5eRYEzihKQStXi1deql0xx32HMHbcoUhOAEAQsK4cXbGa8uW0rff2vFlCCHp6Tblu1Ur6bffpFmzpPff91y3PsEJAOCq3Fzp3nulP/9Z6t/f+poqVXK7KhyRvy3XqJGNav/HP6wZ/Ior3K7MFfQ4AQBcs3u3zUacO1f6978tPCGErFkjDRwozZsn3XyzdeonJLhdlasITgAAVyxfLt1wg53EMWuW3cWOEJGeLg0fLv3rX1L9+tIXX9g+KtiqAwAE33//K3XoIFWpIi1dSmgKGY4jvfOONZj9+9/SM8/Ythyh6QiCE1AEBmACgXP4sN213qePbdHNny/VqeN2VZBkhwF27iz16mWpdt066dFHpZgYtysLKQzABPzEAEzg7KSmSrfcIi1caDtAf/kLQy1DQkZGwbbc734njR0rdevmdlUhix4nAECxW7zYeosPH5a+/lrq2NHtiiDHsXECf/2rtG+f9PTTNm6AFaZTYqsOAFBsHMcWMjp1spuxli0jNIWEH36wo1Fuu01q186GWA4dSmjyA8EJAFAs0tJslWnwYJvTNHeuVKuW21V5XEaGTflu0ULavl367DNp8mQazU4DwQme8Oqrr+qCCy5Q+fLlVb58ebVv314zZ850uywgYn3/vdSmjY0ZmDxZeuklqVQpt6vysPxtucaNpeRk62lavZozbc4AwQmeULt2bY0cOVJLly7V0qVLdfnll+vaa6/V2rVr3S4NiCiOI73+uu3+nHuubc3dcIPbVXncunU25fvWW6WLLrK3H3uMbbkzxF118KxKlSrp+eef1x//+Ee/Pp+76oBTO3BAGjBAevNN6e67rbepTBm3q/KwzExr+P7nP6V69aQxY6QePdyuKuxxVx08x+fz6YMPPlBWVpbat29/0s/Lzs5Wdnb2kbfT09ODUR4QljZskG66Sdq82YJTnz5uV+RhjiN98IHdIbd3rzRsmN05FxvrdmURga06eMbq1at17rnnKiYmRv3799fUqVPVtGnTk37+iBEjFBcXd+SR4PHzmYCTefdd62c6fFj67jtCk6vWr7cp33/4g5SUZHfPPf44oSmA2KqDZ+Tk5Gjr1q3av3+/Jk+erNdee01z5849aXgqbMUpISGBrTrgfzIz7W65iROl22+3EzrKlXO7Ko/KzLTjUV56ye6QGzNGuvJKt6uKSAQneNYVV1yh+vXr6z//+Y9fn0+PE1Bg+XLrNd6xw27S6tOHKeCucBzpww9tW273bmv6fughVpiKEVt18CzHcY5ZUQJQtLw8W9Ro185Wl5Yvl/r2JTS5Yv16Oxrllluk1q1tW+6JJwhNxYzmcHjCY489ph49eighIUEZGRl67733NGfOHH322WdulwaEjdRUqV8/m5n4179K//gHd7S7IivLtuVefNHGsX/8sXTVVW5X5RkEJ3hCamqq7rjjDu3cuVNxcXG64IIL9Nlnn6lr165ulwaEhS++sO04x5FmzpR+/3u3K/Igx7FpokOGSLt2SX/7m/Tww6wwBRk9ToCf6HGCF2Vn2+7P88/brtCbb0rx8W5X5UEbNlgn/qxZUs+e0ujRUmKi21V5Ej1OAIBCrVljg6ZHj5aee85WmghNQZaVZQ3f558vbdokffSRNH06oclFBCcAwDHyG8Bbt5Zyc20200MPSSX4GyN4HEeaMkVq0sT+MB5/XFq7Vrr6arcr8zz+NQCKkJycrKZNmyopKcntUoBit3Wr1KWLNX8PGiQtXSq1bOl2VR7z4492NMqNN0otWtjdck89xfk1IYIeJ8BP9DghkjmO9Pbb0sCBUlycDbW8/HK3q/KYrCzp2WelF16QatWyw/6uucbtqnAcVpwAwON277YTOu64w/qOV60iNAWV40hTp0pNm9qIgaFDbVuO0BSSGEcAAB42darUv7+dM/d//yfdfLPbFXnMxo3SfffZcKwrr5Rmz5bq13e7KpwCK04A4EF79tj5cjfcILVvb200hKYgOnDA5jA1b24TwKdPt0GWhKaQx4oTAHjMtGm2ypSTI02aZAGKI1OCxHEsJA0eLKWkSI8+ag8av8MGK04A4BF790q9e0vXX2/zmdaulXr1IjQFzaZNdjTK9ddbP9OaNdLw4YSmMENwAgAPmDpVatZM+vRT6b//tVWnGjXcrsojDhyQnnzS/gB++MFe/E8+kRo0cLsynAGCEwBEsF9/tQWOG26Q2ra1RY7evVllCor8bbmmTaVRo+xcuR9+kK69lj+AMEZwAorAAEyEo7w86dVX7e/sxYulDz+0hY6aNd2uzCM2b7Yp39ddZ9O/16yR/v53qWxZtyvDWWIAJuAnBmAiXPzwg3TPPdI339jzqFFShQpuV+URBw9KI0fai169ug2xZIUporDiBAARIjtbGjbMjkjZtUuaO1f6z38ITUHz0Ue2xDdypPTgg9K6dbbiRGiKKIwjAIAIsGCBrS5t3Gh3tz/+uBQb63ZVHrF5s3T//dbw3b279MUXUsOGbleFYsKKEwCEsbQ06S9/kTp1sjPmVqywVhpCUxAcPGhLfM2aSatXS1OmSDNnEpoiHCtOABCGHMeOSHngASkjQxo71gJUdLTblXnExx/bUSnbt0sPPSQ99ph0zjluV4UgYMUJAMLM+vVS167SrbdK7dpZM/igQYSmoPjpJzsJ+ZprbGVpzRrpH/8gNHkIK04AECaysqRnnpFefFGqU8eGWfbo4XZVHjBjhjRrlrR7t00SrVZNmjzZBmTR+O05BCcACHGOYzOYBg+WUlPtbNiHH6aPKSj++1+pT5+Ct2+6SZo4kRUmDyM4AUVITk5WcnKyfD6f26XAgzZvlu6913qOr7pK+vpr6Xe/c7uqCHf4sC3nTZxok7/zRUdLCQmEJo9jACbgJwZgIpgyMqQRI6SXXpLi422OYs+e7AwVq1WrLCxNmmSDsC68UGrdWho/3kKTz2dBqmdPtyuFi1hxAoAQkpcnvfWWNHSotH+/9Mgjti3HIkcx2b1beucdC0wrVlj/0h13SH37ShdcYJ9z9dXSnDnSZZcRmsCKE+AvVpxQ3L75xvqYli61O+ZGjbImcATY4cPSZ59JEybYWPCxO0kAAB77SURBVAHHsbvk+vWzbvtSpdyuECGMFScAcNnWrbay9N57tjM0f77UsaPbVUWg1asLtuJ++01q1Up64QXpttukqlXdrg5hguAEAC45cEB67jl7xMVJb7xhO0QlmLAXOLt3S+++a4Fp+XILSL172wvdooXb1SEMEZwAIMgcx/4uf+QRW/gYMsQGT5cr53ZlEeLwYenzzy0szZhhL/jVV0tPPmlbcaVLu10hwhjBCQCCaP58a/ZevNjmJ77wAuMFAmbNmoKtuNRUqWVL6fnnpdtvZysOAUNwAoAgWLfO7pSbPt3ucp89W+rc2e2qIsCePQVbccuWSVWqSL16WaN3y5ZuV4cIRHACgGK0c6c0bJj02mt2h9w770h/+AN9TGclN/fYrbi8PJsO+re/SVdeyVYcihXBCSgCk8NxJjIybBvuhRfsaJQXXpAGDJBiYtyuLIytXWth6b//ta24Cy6wmQ23327zl4AgYI4T4CfmOMEfOTm2ujR8uJSWZnOZHn1UqlDB7crC1N69Nqdh4kRpyRKpcmXbirvzTrbi4ApWnAAgAHJzpbfftm25X36xc2GffpoBlmckN1f64ouCs+J8PtuKmzLFntmKg4sITgBwFvLypA8/lJ56Slq/XrrhBhtG3ayZ25WFoR9+KNiKS0mRzj9fGjnStuKqV3e7OkCSRHsiPGHEiBFKSkpSuXLlVK1aNV133XXasGGD22UhjDmOBaQLL7Rm78REOypl8mRC02nZt0969VXpoovshXv9denmm+0Oue+/lx54gNCEkEJwgifMnTtXAwcO1OLFizVr1izl5uaqW7duysrKcrs0hKGvvpIuvtiON6tQwWYzffqpHZcCP+TmSjNnWuKMj5fuvdfC0eTJ0o4d0pgxlkijotyuFDgBzeHwpF27dqlatWqaO3euLrnkEr++huZwLFokPf649PXXUtu20jPPSFdcwd/vflu3rmArbudOW2G68047AoVVJYQJepzgSWlpaZKkSpUquVwJwsHChdbo/fnndgf89Om22kRg8sO+fdL771tg+vZbqVIl61nq149VJYQlghM8x3EcDRkyRB07dlTz5s1P+nnZ2dnKzs4+8nZ6enowykMImTvXAtPs2bY48t571n7D8Moi+HzSrFkWlqZNs625Hj2si/7qqxlmhbBGcILnDBo0SKtWrdKCBQtO+XkjRozQ8OHDg1QVQoXjWA/T3/8uzZsntWhhrTfXXUdgKtL69QVbcTt2SE2b2n5mr15SjRpuVwcEBD1O8JR7771X06ZN07x585SYmHjKzy1sxSkhIYEepwjlONJnn9kK0+LFUps20pNP2gIJu0mnsH9/wVbc4sVSxYoFW3GtW/PiIeKw4gRPcBxH9957r6ZOnao5c+YUGZokKSYmRjFsKUS8vDwbK/D3v9s4gXbt7Iav7t35O/+kfD7pyy8tLE2dKh0+bFtxH3xgzV/8e4MIRnCCJwwcOFDvvPOOpk+frnLlyiklJUWSFBcXpzJlyrhcHdxw+LD1LI0aZUegdepkbTlduhCYTmrDBunNN6W33pJ+/VVq0sQSZ+/ebMXBM9iqgydEneRvwgkTJqhfv35+XYNxBJEhK8vOknvpJWnrVunKK6VHHpH8nErhPWlpBVtxixbZ4Kr8rbg2bUiZ8BxWnOAJ/P8Bdu+WXn5ZGjvWssBtt0kPP2yneuA4Pp91yOdvxeXk2N7l++9LPXtKsbFuVwi4huAEIKL98ov04ot2kofjSH/6kzRkiFSvntuVhaAffyzYitu+XWrcWBo+3LbiatZ0uzogJBCcAESk776TRo+W/u//pLg46aGHpEGDpCpV3K4sxKSl2Ys0caJN+qxQwZbj+vWTkpLYigOOQ3ACEDFyc21nafRoywC/+531Mv3xj9I557hdXQjx+ezcmAkTpClTbCuuWzfrlr/2WrbigFMgOAFFSE5OVnJysnw+n9ul4CT277eG77FjreH70kttYPXVV0vR0W5XF0I2bizYitu2TWrUSBo2zLbiatVyuzogLHBXHeAn7qoLPZs2SWPGSG+8YYsmt90mDR4stWrldmUhJD29YCvum29s3/LWW+1w3bZt2YoDThMrTgDCiuNIc+ZI//ynDa6sXNmavQcMkOLj3a4uROTl2VbcxIl2XsyhQ7YV9+67thXH7DLgjBGcAISFAwfs7/2xY6Xvv7dDd8ePt5FC5ID/2bTJtuLefNO24s47z86N6d1bql3b7eqAiEBwAhDSNm2SXn3V+pj377eBlS+8wITvIzIyCrbiFiyQype3rbh+/ez8GF4kIKAITgBCjs8nffqp9MordvBupUo2f6l/f7tTzvPy8my/Mn8r7uBBqWtX6Z13pOuuYwkOKEYEJwAhY+dOW1kaP176+Wc70WPCBOkPfyALSJI2by7Yitu6VWrYUHr8cemOO6SEBLerAzyB4ATAVT6f9PnnFpY++kgqXVq65RY73aNtW7erCwEZGdKHH9rq0rx5thX3hz/YVlz79mzFAUFGcALgim3bbIzA66/bP7doYaMFbr/dhld7Wl6eNHeuhaUPP7StuCuukN5+27biypZ1u0LAswhOQBEYgBk4ubnSJ5/Y6tLMmbb9dvvt0t1327ac5xdPfvqpYCvul1+kBg2kxx6zrbg6ddyuDoAYgAn4jQGYZ27LFltZeuMN62NKSrKwdOutUrlyblfnsszMgq24uXPtBcnfirv4YtIkEGJYcQJQLA4csHPjJk6UvvrK8kDv3haYWrZ0uzoXzZghzZ4tVa1qR6B8+KG9WJdfLk2aJF1/PVtxQAgjOAEIGMexUz0mTrTRQhkZ0iWX2ErTzTd7/KDd9HQ7cXj48IL3xcdLjz5qW3F167pXGwC/EZwAnLWtW+3c2DfftIGV9erZMSh9+nh47tKvv9pAyvzHqlXW9J2vRAnbq/zb39yrEcBpIzgBOCNpadKUKba79PXXtrt0003W+H3JJZYLPCMvT1q/3gLS/Pn2/PPP9rEGDaSOHaVBg+wk4gEDpOhom8PQubOrZQM4fQQnAH7Lzra74d5+22Yu5eRIl11mW3E33SSde67bFQZJdra0bFnBatI330h791ogatXKRgZ07Ch16HDiycO1atnU78suk3r2dKN6AGeBu+oAP3n1rrq8PJu7+Pbb1se8f781d/fqZTtNnjg7dv9+aeHCgqD03XcWns45x4ZQduxoj4su8lB6BLyJFScAJ3Aca8l5+23p3Xel7dutb2nAAAtMTZu6XWEx27bt2P6k1avtRale3QLSyJH23KKFVKqU29UCCCKCE1AELw3AXLtW+uADuyNu3TqpcmU7/qRXrwgeKZSXZz/40UFp61b7WKNGFpAeeMCe69eP0BcBgL/YqgP8FKlbdevWWVD6v/+TfvjBjkK79loLTN262dlxEeXQIWnp0oJG7oULbSuuZEnpwgsLtt06dJCqVXO7WgAhhhUnwIM2bCgIS2vW2HDKa6+1Hahu3aSYGLcrDKC9e4/tT1qyxLrazz3XltGGDLGg1LatxwdNAfAHK06An8J5xclxbGVpyhTbilu1ynJDz562stS9uxQb63aVAeA4dsbb0dtua9fax2rUsIDUqZM9n3++rTIBwGkgOAF+CrfglJdniytTp9rjxx9tQeXosFSmjNtVniWfz5bMjg5K27fbx5o0Kdh269hRSkykPwnAWeN/t4AIcviwnRM7dao0bZq0Y4dUpYqFpRdflK64IsxXlg4etFEA+SFp4UI7yqRUKalNG+m22ywkXXyx/eAAEGAEJyDMHTggff65haWPP5b27ZPq1LGz4a6/3nqcw3ZHavduC0f507iXLbN0WL68haOHH7aglJTEwbgAgiJc/3MKeNq2bdInn9jjq69sIaZZM2ngQAtLrVqF4a6U40hbthy77bZunX2sVi3rTerd24JS8+Y2pRsAgoweJ8BPbvY4+Xy2Q/XxxxaWvv/eckPHjtJVV9kJHw0bBrWks5eba13qRwelnTvtY82aFfQmdepkS2hhlwQBRCJWnIAiuDUAMy3NtuA++UT69FPbtapcWerRQxo61Jq7K1QIaklnJyvr2P6kRYukjAwbFJWUJPXpU9CfVKmS29UCQKFYcQL8VNwrTo5jd77lryrNn2+LMhdcYKtKV19tR6GFzQ7Vb7/Z4bf5QWn5cvuBKlSwxqv8FaU2bcK8Yx2Al7DiBLgoLU2aPdtWlr74wlp8YmOlyy+XxoyxwFSnjttV+sFxpM2bC6ZxL1hgKVCyH6BjR6lfP3tu1kwqUcLVcgHgTBGcgCDy+ey0jy++sLC0eLG9r2FDC0ndu1toCvkbxHJzpZUrj+1PSk21PqTzz7e5B8OG2cpSWCQ/APAPwQkoZtu2FQSlL7+0cQHly0tdukjJyXbESWKi21UWITNT+vbbgtWkxYutZykmxo4quesua+Ju3z7MGq8A4PTQ4wTPmDdvnp5//nktW7ZMO3fu1NSpU3Xdddf5/fX+9jjt3y/Nm2dbcF98YXfUlyhh/c/du1tQuuiiEJ+tlJJybH/SihW2NFax4rHTuFu3jrCD7QDg1EL5P91AQGVlZalFixa68847deONNwbwupYtZs+2x/LldtxJnTpS167S8OG2uhSyN4rld6Ufve22aZN9LDHRAtLdd9tz48b0JwHwNIITPKNHjx7q0aPHWV8nO9uONckPSt9+a8Os4+OtP6l/f6lz5xA+Gu3wYVtBOjoo7dplxbZoIf3+9xaSOnSQatd2u1oACCkEJ+AksrOzlZ2drZwcaeXKaM2adUiSlJBg4alSJQtI//ynBabGjUM0KKWnW09SfkhavNhGjcfGSu3aSX/+swWldu2kuDi3qwWAkEZwAo6TmWmzGZ95ZrHmzcuTdJGkspL2SZKefFK68kqbrxSSu1Y7dlh/Un4j9/ff295h5coWkJ5+2hq5W7Wy4ZMAAL8RnOB5u3YVjB+aP7+gD7pKlUt09dW5at/ep/btM5WYaNtvgwbZXXGumzHD9gobNrQAlL+i9NNP9vH69S0oDRhgz40aheiSGACED+6qg6c4jvTLL1Ji4h3q2vXv2ratntavt4/VrWsLMfmP47fe3DyrTjk5Nh3zxx/tMWuWzTfIFxVlK0idOhX0J9WoEdwaAcADCE6IaNafZG09ixbZgsz27faxOnXSdOWVcUeCUkLCqa9V7MHJ55O2bpU2brRwdPTzli223SbZdMyyZaU9eywJligh/eUv0ssvB74mAMAx2KpDRNmxwwJSflBatkw6dMhGDbVo4VOXLrvVqlWWBg9uo8GDn1Dnzp1VqVIlJSQEabq140g7dxYejjZtsqQnSaVK2VbbeedJ111nzw0b2nPNmtJHH0nXXmsH1/l8NhwKAFDsWHFC2MrOttWkRYsKwtLWrfaxOnVsiHW7dvbcsqW0aNEcde7c+YTr9O3bVxMnTizy+53WitOePSeGo/x/zsqyzylRwvYHjw5F+c916hQ9IXPGDGnOHOmyy6SePYusHwBw9ghOCAu5udIPP0hLlthZb0uWSKtW2UiimBipTZtjg1LNmoGv4YTglJl5YijKf967t+ALa9YsPBz97ndM3QaAMMNWHUJOXp5lj/yAtGSJ3el28KD1QDdtakHpzjvtGJOWLYvxrvpDh6TNm62gVavsfT16WM/Rzp0Fn1e5soWhRo2ka64pCEgNGkjnnltMxQEAgo3gBFc5jm2v5QekpUutLyktzT5ev76FoxtvtLB04YXFkENyc6Wffy687+iXX6xIqeAb16hhZ6jkh6OGDUP4PBUAQCCxVYeg8fmkDRts9WjlSntesaJgV6t2bQtHSUn2aN06gHkkL0/69dfC+45++snCk2RbZw0aHLOlNmXNGo3+5BOlSvpx40Z3xhEAAEICwQnF4uBBafXqY0PSqlX2fkmqV8+22Fq1skdSkp31dlYcx6ZZFtZ3tGlTwTePjrZJluedd2LvUULCSceBuzrHCQAQEghOOCv5d9evXm2P77+3kLR+va0wRUfbIMn8gNSqlQWmihXP4pumpZ24pZb/nL/HJ1kIKiwcJSba7f6nieAEACA4wW/p6dKaNfbID0qrVxdstZUtK51//rEhqXlzqUyZM/hmBw7YKlFh4ei33wo+r3r1E+9Wa9jQmqPKlg3Iz52P4AQAIDjhBIcPWy/S0eFo9Wrrk5ZsFalhQwtJRz8SE0/z0Nv8Y0QK6zvKH+8tSXFxBStHxwekIAYYghMAgODkYYcPW1/0unX2yF9JWr/ePiZJtWqdGJAaN5ZiY/38Jj6ftG1b4eHo55/t45ItSxU266hhQ6lKlZA4nJbgBAAgOHlAVpatIOUHpPzHpk0FASkuzrbVzj//2Ge/7mpzHCklpfCm7M2bbcS3VHCMSGHhqGbN01yuCj6CEwCA4BRBdu2y1aLjA1L+MSSSrSA1aWKrRk2aFDyqV/djUWfv3pNPys7MtM+JirJb5o4PR/4eIxLCCE4AAIJTmDl4sKBn+ugMs26dHY8mWQ9S/fonBqTGjU/REjRjhvT113ZmyXnnFd6Unf8NJE8eI0JwAgCE7//+R7CT9Uxv3GjtQvnKly/ILN26FQSkBg0KyS45OXY32sZUKTXVttbyn5cvlxYsOLGQ/GNEzjtPuvrqgnDksWNEkpOTlZycLF9+PxYAwLNYcXLJwYMWjjZvtgbt/OPQNm48sWc6/1SP4xd4qlbMVdTuXceGoJM9H33obL7KlW3qZHq63cXmONZndMstUnIyx4gchxUnAAArTsXEcWyBJz8Y5Yej/Oejz4eNibFb+Rs2lK7v6dP5NXarccUUJZZJVeXDKSqx63+rRFtTpCVHBaLduwvOUctXsaI1LMXH2/MFFxz7dv5ztWoFQyBnzJCuvdb2+Hw+6bbbCE0AABSCFaezkJ5us41++cVWiY4PSAcO2OdFKU+Nq+xRq5qpalYlReeVT1XdmBTFR1kwKpOeqqjU/4WhXbvsXLWjxcUVHn6Of1+1amfeXzRjhjRnjnTZZVLPnmfxqkQuVpwAAASnk3AcW9DJD0YnPH52FLV/r6orVfFKUa3oVDWplKIG56aqTukUVVeqKuak6JzMVJXc+5uiju+POffcE4NQYc/Vqp3h6G0EGsEJAODZrTqfz7bLjg9Duzen6cBPKfLtSFWF7JRjglHr2BTViE5VVV+Kyh/6TdE6fNQFJWWVlcrFS5Xyw0+7wleJqleXzjnHtZ8dAACcmYgMTrm51ga0fbv063ZHv23OUPqPKTr4c6p8O1JU4rdUxexLUVXHQlEDpapTiVRVd1JU2sk55lp5MbFyqlVXiZrxiqpeXYpvffLVIQ/daQYAgBeFXXDKzpZ27JB2/JipvT+kKGNTqg79kirfrymK2pWq2H0pKncgVdVlq0UtlKIyOnTMNXKjS+tgheryVYlXiZrVFVuvhUrXLryHqET58iFx3AcAAHBfyASn7Gwp5acD2vNDqtJ/TNGBLak6vM0apkvuTlFseqrKZ6WoUq6tEiXqwDFfnxtVUullqutQ+erK/V28omo0VUydy5VXv7qcxHhFxRcEopIVKqgcYQgAAJymwAen/AnUnTtLPXsqY9ch7VqTqv3rU5S5OVXZW1OVt8NusY/ZZ83TFbJTVDUvVXWVobpHXSpX0dpXqpoyylTXgbh4+RLPU1aNS7S9TnWd2yBeFRpV1zn1/xeGKlZUpRA/6wwAAIS3wN5V9795QI6kKEmZKqtzj1sZylOU9kRXU1pMdWWVq67sivFyqloPUcxRgahC43iVqFo55A9+ReQ7enL4jz/+yF11AOBhgQ1ODzygvNH/Ugk5ylOUdiR21Lbuf1LZ38Ur7rzqqtwsXufWq6KoktEB+5ZAsDCOAAAQ2K26zp1VYvRoKTpaJXw+1R79oGozTBEAAESIwAannj2l6dOZQA0AACISk8MBP7FVBwCg8xoAAMBPBCcAAAA/EZwAAAD8RI8T4CfHcZSRkaFy5copisnzAOBJBCcAAAA/sVUHAADgJ4ITAACAnwhOAAAAfiI4AQAA+IngBAAA4CeCEwAAgJ8ITgAAAH4iOAEAAPiJ4AQAAOAnghMAAICfSvrzSflndAEAAEQqf84i9Ss4ZWRkKC4uLiBFAQAAhKK0tDSVL1/+lJ/j1yG/p7PilJ6eroSEBG3btq3Ib16UpKQkLVmy5KyuEehrhdJ1AvlaB6qmSL5OKP5uh9p1AnWtUHytA3mtULoO/x0J7nVC8Xc71K4TqGudyWsdsBWnqKio0/4DLl++/Fn/UkRHRwfkX+RAXivUriMF5rWWQu9nC7Xr5Aul3+1Qu06grxVKr3UgrxVq15H470iwrpMvlH63Q+06gb5WoH6384V0c/jAgQND7lqhdp1ACrWfLdSuE0ih9rOF4r9rgRKKP1uoXSeQQu1nC7XrBFKo/Wyh+O9acfBrq+50pKenKy4uzq99QpwdXuvg4vUOHl7r4OG1Di5e7+Aprtc6etiwYcMCdrX8i0ZH67LLLlPJkn7tBOIs8FoHF6938PBaBw+vdXDxegdPcbzWAV9xAgAAiFQh3eMEAAAQSghOAAAAfiI4AQAA+IngBAAA4KeAB6cpU6aoe/fuqlKliqKiorRy5cpAfwtPcRxHw4YNU82aNVWmTBlddtllWrt27Sm/ZuLEiYqKijrhcejQoSBVHXleeeUVJSYmKjY2Vq1bt9b8+fPdLinsnc5rOmfOnEJ/p9evXx/EiiPLvHnzdM0116hmzZqKiorStGnT3C4p7J3ua8rvdfEYMWKEkpKSVK5cOVWrVk3XXXedNmzYELDrBzw4ZWVlqUOHDho5cmSgL+1Jzz33nF566SW9/PLLWrJkieLj49W1a9cij8ApX768du7cecwjNjY2SFVHlvfff1+DBw/W448/rhUrVqhTp07q0aOHtm7d6nZpYetMX9MNGzYc8zvdsGHDIFUcebKystSiRQu9/PLLbpcSMc70NeX3OrDmzp2rgQMHavHixZo1a5Zyc3PVrVs3ZWVlBeYbOMVky5YtjiRnxYoVxfUtIl5eXp4THx/vjBw58sj7Dh065MTFxTn//ve/T/p1EyZMcOLi4oJRoie0bdvW6d+//zHva9y4sfPoo4+6VFH4O93X9Ouvv3YkOfv27QtGeZ4jyZk6darbZUQUf15Tfq+D47fffnMkOXPnzg3I9ehxCmFbtmxRSkqKunXrduR9MTExuvTSS7Vw4cJTfm1mZqbq1q2r2rVr6+qrr9aKFSuKu9yIlJOTo2XLlh3zZyBJ3bp1K/LPAIU7m9e0VatWqlGjhrp06aKvv/66OMsEgobf6+KVlpYmSapUqVJArkdwCmEpKSmSpOrVqx/z/urVqx/5WGEaN26siRMnasaMGXr33XcVGxurDh06aOPGjcVabyTavXu3fD7faf8Z4OTO5DWtUaOGxo0bp8mTJ2vKlClq1KiRunTponnz5gWjZKBY8Htd/BzH0ZAhQ9SxY0c1b948INc8qxnkb7/9tv785z8feXvmzJnq1KnTWRflVce/np988okkKSoq6pjPcxznhPcdrV27dmrXrt2Rtzt06KALL7xQY8eO1ZgxYwJctTec7p8BinY6r2mjRo3UqFGjI2+3b99e27Zt0wsvvKBLLrmkWOsEigu/18Vv0KBBWrVqlRYsWBCwa55VcOrZs6cuuuiiI2/XqlXrrAvysuNfz+zsbEm28lSjRo0j7//tt99O+L/1UylRooSSkpJYcToDVapUUXR09AkrIaf7Z4ACgXpN27Vrp0mTJgW6PMBV/F4Hzr333qsZM2Zo3rx5ql27dsCue1ZbdeXKlVODBg2OPMqUKROoujzp+NezadOmio+P16xZs458Tk5OjubOnauLL77Y7+s6jqOVK1ceE77gn9KlS6t169bH/BlI0qxZs07rzwAFAvWarlixgt9pRBx+r8+e4zgaNGiQpkyZotmzZysxMTGg1w/40cx79+7V1q1btWPHDkk6MjshPj5e8fHxgf52ES0qKkqDBw/Ws88+q4YNG6phw4Z69tlnVbZsWd1+++1HPq9Pnz6qVauWRowYIUkaPny42rVrp4YNGyo9PV1jxozRypUrlZyc7NaPEtaGDBmiO+64Q23atFH79u01btw4bd26Vf3793e7tLBV1Gs6dOhQ/frrr3rrrbckSaNHj1a9evXUrFkz5eTkaNKkSZo8ebImT57s5o8R1jIzM7Vp06Yjb2/ZskUrV65UpUqVVKdOHRcrC19Fvab8XgfHwIED9c4772j69OkqV67ckdXtuLi4wCzwBOTevKNMmDDBkXTC46mnngr0t/KEvLw856mnnnLi4+OdmJgY55JLLnFWr159zOdceumlTt++fY+8PXjwYKdOnTpO6dKlnapVqzrdunVzFi5cGOTKI0tycrJTt25dp3Tp0s6FF14YsNtavexUr2nfvn2dSy+99Mjbo0aNcurXr+/ExsY6FStWdDp27Oh88sknLlQdOfJvhT/+cfR/S3B6inpN+b0OjsL+DCQ5EyZMCMj1o/73TQAAAFAExhEAAAD4ieAEAADgJ4ITAACAnwhOAAAAfiI4AQAA+IngBAAA4CeCEwAAgJ8ITgAAAH4iOAEAAPiJ4AQAAOAnghMAAICfCE4AAAB++n9rbrq56HnpOgAAAABJRU5ErkJggg==\n",
148 "text/plain": [
149 "Graphics object consisting of 2 graphics primitives"
150 ]
151 },
152 "execution_count": 4,
153 "metadata": {},
154 "output_type": "execute_result"
155 }
156 ],
157 "source": [
158 "var('y')\n",
159 "\n",
160 "def euler_desolve(f, x0, y0, x1):\n",
161 " n = 5\n",
162 " h = (x1-x0)/n\n",
163 " S = []\n",
164 " Y = [y0]\n",
165 " for i in range(n+1):\n",
166 " S.append(x0 + i*h)\n",
167 " Y.append(N( Y[i] + h*f(S[i], Y[i]) ))\n",
168 " return S, Y\n",
169 "\n",
170 "f(x,y) = y\n",
171 "x0 = -1\n",
172 "x1 = 2\n",
173 "y0 = e^(-1)\n",
174 "\n",
175 "S, Y = euler_desolve(f, x0, y0, x1)\n",
176 "plot(e^x, -1, 2) + line([(S[i], Y[i]) for i in range(len(S))], color='red', marker='o', markersize=2)"
177 ]
178 },
179 {
180 "cell_type": "markdown",
181 "metadata": {},
182 "source": [
183 "Sage also has an `eulers_method()` function \"for pedagogical purposes only\":"
184 ]
185 },
186 {
187 "cell_type": "code",
188 "execution_count": 5,
189 "metadata": {},
190 "outputs": [
191 {
192 "name": "stdout",
193 "output_type": "stream",
194 "text": [
195 " x y h*f(x,y)\n",
196 " -1 0.367879441171442 0.0367879441171442\n",
197 "-0.900000000000000 0.404667385288587 0.0404667385288587\n",
198 "-0.800000000000000 0.445134123817445 0.0445134123817445\n",
199 "-0.700000000000000 0.489647536199190 0.0489647536199190\n",
200 "-0.600000000000000 0.538612289819109 0.0538612289819109\n",
201 "-0.500000000000000 0.592473518801020 0.0592473518801020\n",
202 "-0.400000000000000 0.651720870681122 0.0651720870681122\n",
203 "-0.300000000000000 0.716892957749234 0.0716892957749234\n",
204 "-0.200000000000000 0.788582253524157 0.0788582253524157\n",
205 "-0.100000000000000 0.867440478876573 0.0867440478876573\n",
206 "-1.38777878078145e-16 0.954184526764230 0.0954184526764230\n",
207 "0.0999999999999999 1.04960297944065 0.104960297944065\n",
208 "0.200000000000000 1.15456327738472 0.115456327738472\n",
209 "0.300000000000000 1.27001960512319 0.127001960512319\n",
210 "0.400000000000000 1.39702156563551 0.139702156563551\n",
211 "0.500000000000000 1.53672372219906 0.153672372219906\n",
212 "0.600000000000000 1.69039609441897 0.169039609441897\n",
213 "0.700000000000000 1.85943570386086 0.185943570386086\n",
214 "0.800000000000000 2.04537927424695 0.204537927424695\n",
215 "0.900000000000000 2.24991720167165 0.224991720167165\n",
216 "1.00000000000000 2.47490892183881 0.247490892183881\n",
217 "1.10000000000000 2.72239981402269 0.272239981402269\n",
218 "1.20000000000000 2.99463979542496 0.299463979542496\n",
219 "1.30000000000000 3.29410377496746 0.329410377496746\n",
220 "1.40000000000000 3.62351415246420 0.362351415246420\n",
221 "1.50000000000000 3.98586556771062 0.398586556771062\n",
222 "1.60000000000000 4.38445212448168 0.438445212448168\n",
223 "1.70000000000000 4.82289733692985 0.482289733692985\n",
224 "1.80000000000000 5.30518707062284 0.530518707062284\n",
225 "1.90000000000000 5.83570577768512 0.583570577768512\n",
226 "2.00000000000000 6.41927635545363 0.641927635545363\n"
227 ]
228 }
229 ],
230 "source": [
231 "# Usage: eulers_method(f, x0, y0, h, x1)\n",
232 "eulers_method(f, -1, N(e^(-1)), 0.1, 2)"
233 ]
234 },
235 {
236 "cell_type": "markdown",
237 "metadata": {},
238 "source": [
239 "## Solving the heat equation with a finite difference method"
240 ]
241 },
242 {
243 "cell_type": "code",
244 "execution_count": null,
245 "metadata": {},
246 "outputs": [],
247 "source": [
248 "def heat_fdm(u0j, u1j, ui0):\n",
249 " m, n = len(u0j)-1, len(ui0)-1\n",
250 " k, h = 1/m, 1/n\n",
251 " \n",
252 " u = [[0] * (m+1) for i in range(n+1)]\n",
253 " for j in range(m+1):\n",
254 " u[0][j] = u0j[j]\n",
255 " for j in range(m+1):\n",
256 " u[n][j] = u1j[j]\n",
257 " for i in range(n+1):\n",
258 " u[i][0] = ui0[i]\n",
259 " \n",
260 " for j in range(0,m):\n",
261 " for i in range(1,n):\n",
262 " u[i][j+1] = (k/(h*h)) * (u[i+1][j] - 2*u[i][j] + u[i-1][j]) + u[i][j]\n",
263 " \n",
264 " return u\n",
265 "\n",
266 "n, m = 20, 20\n",
267 "u0j = [10 - (j/m)*10 for j in range(m+1)] # One extreme goes from hot to cold\n",
268 "u1j = [(j/m)*10 for j in range(m+1)] # The other does the opposite\n",
269 "ui0 = [10 - (i/m)*10 for i in range(0,n+1)]\n",
270 "\n",
271 "u = heat_fdm(u0j, u1j, ui0)\n",
272 "for t in range(m+1):\n",
273 " show(line([(i/n, u[i][t]) for i in range(n+1)], ymin=-1, ymax =12))"
274 ]
275 }
276 ],
277 "metadata": {
278 "kernelspec": {
279 "display_name": "SageMath 9.0",
280 "language": "sage",
281 "name": "sagemath"
282 },
283 "language_info": {
284 "codemirror_mode": {
285 "name": "ipython",
286 "version": 3
287 },
288 "file_extension": ".py",
289 "mimetype": "text/x-python",
290 "name": "python",
291 "nbconvert_exporter": "python",
292 "pygments_lexer": "ipython3",
293 "version": "3.8.5"
294 }
295 },
296 "nbformat": 4,
297 "nbformat_minor": 4
298}

Generated with cgit - Back to sebastiano.tronto.net