{ "cells": [ { "cell_type": "markdown", "id": "3715d656", "metadata": {}, "source": [ "# Integralrechnung" ] }, { "cell_type": "markdown", "id": "fe3eefe5", "metadata": {}, "source": [ "*Notebook erstellt am 29.09.2022 von C. Rockstuhl, überarbeitet von Y. Augenstein*\n", "\n", "*Abbildungen in diesem Notebook wurden dem Script \"Computational Physics\" von Prof. T. Pertsch (Uni-Jena) entnommen*" ] }, { "cell_type": "markdown", "id": "953d0964", "metadata": {}, "source": [ "In diesem Notebook werden wir einige Details zur Integralrechnung besprechen. Wir werden dieses Notebook in zwei Abschnitte unterteilen und einmal diskutieren, wie wir symbolisch ausgewählte Integrale lösen können, und zum anderen werden wir diskutieren, wie wir diese numerisch evaluieren können. " ] }, { "cell_type": "markdown", "id": "e645a4fa", "metadata": {}, "source": [ "## Symbolisches Lösen von Integralen " ] }, { "cell_type": "markdown", "id": "43ee9330", "metadata": {}, "source": [ "Wie bereits in früheren Notebooks, in denen wir symbolisch differenziert haben, benutzen wir hier wieder die Bibliothek SymPy. Mehr Informationen und eine ausführliche Dokumentation finden Sie hier: [https://www.sympy.org](https://www.sympy.org/en/index.html).\n", "\n", "Hierfür importieren wir wieder das Modul `sympy` sowie einige symbolische Variablen, die wir noch benötigen werden." ] }, { "cell_type": "code", "execution_count": 1, "id": "6a5a54d3", "metadata": {}, "outputs": [], "source": [ "from IPython.display import display # zum Anzeigen von Gleichungen\n", "import sympy\n", "from sympy.abc import a, b, x, y\n", "from sympy import Integral, Eq, S # dient nur der Anzeige von Integralen und Gleichungen" ] }, { "cell_type": "markdown", "id": "1892783f", "metadata": {}, "source": [ "Nun definieren wir wieder eine einfache Funktion. Wir nehmen hier ein Beispiel aus der Vorlesung." ] }, { "cell_type": "code", "execution_count": 2, "id": "002d24a9", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle f = 3 x^{3} - 2 x$" ], "text/plain": [ "Eq(f, 3*x**3 - 2*x)" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f = 3 * x**3 - 2 * x\n", "\n", "Eq(S(\"f\"), f) # Dient nur der Ausgabe. Die Logik der Ausgabe ist, dass die linke und die rechte Seite der \n", " # Gleichung spezifiziert werden und man sich diese Anzeigen läst. \n", " # Die Funktion S 'sympifiziert' einmal den vorhandenen Ausdruck. In unserem Fall einfach f.\n", " " ] }, { "cell_type": "markdown", "id": "0d0cb6f7", "metadata": {}, "source": [ "Das Erste, was wie wieder machen möchten, ist die Evaluation der Funktion an einem bestimmten Wert für $x$. Hierfür können wie den Befehl `f.subs` verwenden, der einen entsprechenden Wert für $x$ substituiert in den Ausdruck der Funktion $f(x)$." ] }, { "cell_type": "code", "execution_count": 3, "id": "e3002298", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3(2)^3 - 2*2 = 20\n" ] } ], "source": [ "z = f.subs(x, 2)\n", "print(f\"3(2)^3 - 2*2 = {z}\") # Und wir geben uns das Ergebnis aus." ] }, { "cell_type": "markdown", "id": "b7e51b28", "metadata": {}, "source": [ "Um jetzt unsere Funktion `f` bezüglich `x` zu integrieren, verwenden wir `integrate()`." ] }, { "cell_type": "code", "execution_count": 4, "id": "5b4e05df", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int \\left(3 x^{3} - 2 x\\right)\\, dx = \\frac{3 x^{4}}{4} - x^{2}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, x), 3*x**4/4 - x**2)" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = f.integrate(x)\n", "\n", "Eq(Integral(f, x), intf) # Hier dient die Funktion Integral zur Visulisierung des entsprechenden Integrals.\n", " # Und auch hier werden wieder die linke und die rechte Seite der Gleichung \n", " # spezifiziert." ] }, { "cell_type": "markdown", "id": "7b9e49c2", "metadata": {}, "source": [ "Alternativ kann man für diesen Zweck auch den etwas allgemeineren Funktionsaufruf verwenden." ] }, { "cell_type": "code", "execution_count": 5, "id": "cd3df46b", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int \\left(3 x^{3} - 2 x\\right)\\, dx = \\frac{3 x^{4}}{4} - x^{2}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, x), 3*x**4/4 - x**2)" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = sympy.integrate(f, x)\n", "\n", "Eq(Integral(f, x), intf)" ] }, { "cell_type": "markdown", "id": "bc59fd69", "metadata": {}, "source": [ "Die obigen Beispiele sind unbestimmte Integrale von $f(x)$. Wir können aber auch bestimmte Integrale mit derselben `integrate()`-Funktion auswerten. Wir müssen hierfür lediglich der Funktion die entsprechenden Integrationsgrenzen als Argumente mit übergeben." ] }, { "cell_type": "code", "execution_count": 6, "id": "6a3c5813", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int\\limits_{0}^{3} \\left(3 x^{3} - 2 x\\right)\\, dx = \\frac{207}{4}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, (x, 0, 3)), 207/4)" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = sympy.integrate(f, (x, 0, 3)) # alternativ auch hier wieder: f.integrate((x, 0, 3))\n", "\n", "Eq(Integral(f, (x, 0, 3)), intf)" ] }, { "cell_type": "markdown", "id": "c711096f", "metadata": {}, "source": [ "Beachten Sie bitte, dass Sie an dieser Stelle immer noch rein symbolisch rechnen. Wenn Sie das gewünschte Ergebnis in einer numerisch besser verarbeitbaren Form benötigen, müssen Sie es z.B. in eine Gleitkommazahl (`float`) umwandeln. Mehr Informationen zu numerischen Typen finden Sie hier: https://docs.python.org/3/library/stdtypes.html#numeric-types-int-float-complex" ] }, { "cell_type": "code", "execution_count": 7, "id": "28ba293d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "51.75" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "float(intf)" ] }, { "cell_type": "markdown", "id": "dd2588b8", "metadata": {}, "source": [ "Ebenso können Sie Symbole als Integrationsgrenzen verwenden, um z.B. einfach nur zwischen $a$ und $b$ zu integrieren." ] }, { "cell_type": "code", "execution_count": 8, "id": "9a020149", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int\\limits_{a}^{b} \\left(3 x^{3} - 2 x\\right)\\, dx = - \\frac{3 a^{4}}{4} + a^{2} + \\frac{3 b^{4}}{4} - b^{2}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, (x, a, b)), -3*a**4/4 + a**2 + 3*b**4/4 - b**2)" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = f.integrate((x, a, b))\n", "\n", "Eq(Integral(f, (x, a, b)), intf)" ] }, { "cell_type": "markdown", "id": "1a125578", "metadata": {}, "source": [ "Wir können unsere Funktionen auch komplizierter machen und sie, zum Beispiel, von zwei Variablen abhängen lassen. In unserem Fall wäre das $x$ und $y$ und die Funktion lautet:\n", "$$\n", "g(x,y)=x \\sin (xy). \n", "$$" ] }, { "cell_type": "code", "execution_count": 9, "id": "5d6e79e7", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle g = x \\sin{\\left(x y \\right)}$" ], "text/plain": [ "Eq(g, x*sin(x*y))" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g = x * sympy.sin(x * y)\n", "\n", "Eq(S(\"g\"), g)" ] }, { "cell_type": "code", "execution_count": 10, "id": "43b34680", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int x \\sin{\\left(x y \\right)}\\, dx = \\begin{cases} - \\frac{x \\cos{\\left(x y \\right)}}{y} + \\frac{\\sin{\\left(x y \\right)}}{y^{2}} & \\text{for}\\: y \\neq 0 \\\\0 & \\text{otherwise} \\end{cases}$" ], "text/plain": [ "Eq(Integral(x*sin(x*y), x), Piecewise((-x*cos(x*y)/y + sin(x*y)/y**2, Ne(y, 0)), (0, True)))" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intg = g.integrate(x)\n", "\n", "Eq(Integral(g, x), intg)" ] }, { "cell_type": "markdown", "id": "e4995663", "metadata": {}, "source": [ "Beachten Sie, dass das Ergebnis dieses Integrals davon abhängt, ob $y$ Null ist oder eben gerade verschieden von Null. Wir können auch von Anfang an spezifizieren und bei der erstmaligen Definition festlegen, dass unsere Variablen ungleich Null sind." ] }, { "cell_type": "code", "execution_count": 11, "id": "4b5bdc47", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle g = x \\sin{\\left(x y \\right)}$" ], "text/plain": [ "Eq(g, x*sin(x*y))" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x = sympy.Symbol(\"x\", nonzero=True)\n", "y = sympy.Symbol(\"y\", nonzero=True)\n", "\n", "g = x * sympy.sin(x * y)\n", "\n", "Eq(S(\"g\"), g)" ] }, { "cell_type": "markdown", "id": "e94febf3", "metadata": {}, "source": [ "Wenn wir nun wieder wie oben über $x$ integrieren, erhalten wir nur noch ein Ergebnis." ] }, { "cell_type": "code", "execution_count": 12, "id": "3aad1d07", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int x \\sin{\\left(x y \\right)}\\, dx = - \\frac{x \\cos{\\left(x y \\right)}}{y} + \\frac{\\sin{\\left(x y \\right)}}{y^{2}}$" ], "text/plain": [ "Eq(Integral(x*sin(x*y), x), -x*cos(x*y)/y + sin(x*y)/y**2)" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intgx = g.integrate(x)\n", "\n", "Eq(Integral(g, x), intgx)" ] }, { "cell_type": "markdown", "id": "1b28173a", "metadata": {}, "source": [ "Wir können selbstverständlich auch über die andere Koordinate ($y$) integrieren." ] }, { "cell_type": "code", "execution_count": 13, "id": "6feaf6df", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int x \\sin{\\left(x y \\right)}\\, dy = - \\cos{\\left(x y \\right)}$" ], "text/plain": [ "Eq(Integral(x*sin(x*y), y), -cos(x*y))" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intgy = g.integrate(y)\n", "\n", "Eq(Integral(g, y), intgy)" ] }, { "cell_type": "markdown", "id": "786f11ff", "metadata": {}, "source": [ "Um einen etwas allgemeineren Zugang zu bekommen, können wir im Folgenden auch die eigentliche Funktionsdefinition in eine Python Funktion auslagern und dann im Folgenden nur mit einem allgemeinen Beispiel rechnen. Dafür definieren wir uns als erstes eine passende Funktion. " ] }, { "cell_type": "code", "execution_count": 14, "id": "2071f226", "metadata": {}, "outputs": [], "source": [ "def f(x):\n", " return 3 * x**3 - 2 * x" ] }, { "cell_type": "markdown", "id": "70057b97", "metadata": {}, "source": [ "Und dann können wir alle Befehle, die wir oben verwendet haben, ganz allgemein auf die Funktion $f$ anwenden und ganz normal damit rechnen. " ] }, { "cell_type": "code", "execution_count": 15, "id": "afa45dfe", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int \\left(3 x^{3} - 2 x\\right)\\, dx = \\frac{3 x^{4}}{4} - x^{2}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, x), 3*x**4/4 - x**2)" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = f(x).integrate(x)\n", "\n", "Eq(Integral(f(x), x), intf)" ] }, { "cell_type": "code", "execution_count": 16, "id": "1064483a", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int\\limits_{a}^{b} \\left(3 x^{3} - 2 x\\right)\\, dx = - \\frac{3 a^{4}}{4} + a^{2} + \\frac{3 b^{4}}{4} - b^{2}$" ], "text/plain": [ "Eq(Integral(3*x**3 - 2*x, (x, a, b)), -3*a**4/4 + a**2 + 3*b**4/4 - b**2)" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = sympy.integrate(f(x), (x, a, b))\n", "\n", "Eq(Integral(f(x), (x, a, b)), intf)" ] }, { "cell_type": "markdown", "id": "4fee4b9a", "metadata": {}, "source": [ "Wir können genauso gut Funktionen mit mehreren Variablen verwenden und dann anschließend doppelte Integrale durchführen.\n", "\n", "Hierfür definieren wir zunächst eine Funktion, die von $x$ und $y$ abhängt." ] }, { "cell_type": "code", "execution_count": 17, "id": "e3329e05", "metadata": {}, "outputs": [], "source": [ "def g(x, y):\n", " return x * sympy.sin(x * y)" ] }, { "cell_type": "markdown", "id": "9bdd5ead", "metadata": {}, "source": [ "Das unbestimmt Integral können wir nun auswerten, indem wir explizit zwei Integrale ausführen." ] }, { "cell_type": "code", "execution_count": 18, "id": "02a3d3f9", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\iint x \\sin{\\left(x y \\right)}\\, dy\\, dx = - \\frac{\\sin{\\left(x y \\right)}}{y}$" ], "text/plain": [ "Eq(Integral(x*sin(x*y), y, x), -sin(x*y)/y)" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intg = sympy.integrate(sympy.integrate(g(x, y), y), x)\n", "\n", "Eq(Integral(Integral(g(x, y), y), x), intg)" ] }, { "cell_type": "markdown", "id": "52cb3578", "metadata": {}, "source": [ "Oder alternativ können wir auch in einem Aufruf von `integrate()` mehrere Variablen angeben." ] }, { "cell_type": "code", "execution_count": 19, "id": "82ee9791", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\iint x \\sin{\\left(x y \\right)}\\, dy\\, dx = - \\frac{\\sin{\\left(x y \\right)}}{y}$" ], "text/plain": [ "Eq(Integral(x*sin(x*y), y, x), -sin(x*y)/y)" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intg = sympy.integrate(g(x, y), y, x)\n", "\n", "Eq(Integral(g(x, y), y, x), intg)" ] }, { "cell_type": "markdown", "id": "08730213", "metadata": {}, "source": [ "Bestimmte Integrale können wir ganz analog auswerten." ] }, { "cell_type": "code", "execution_count": 20, "id": "60745a3e", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int\\limits_{-2}^{-1}\\int\\limits_{1}^{3} x \\sin{\\left(x y \\right)}\\, dy\\, dx = - \\sin{\\left(1 \\right)} + \\frac{\\sin{\\left(3 \\right)}}{3} - \\frac{\\sin{\\left(6 \\right)}}{3} + \\sin{\\left(2 \\right)} = 0.20800494410405$" ], "text/plain": [ "Eq(Eq(Integral(x*sin(x*y), (y, 1, 3), (x, -2, -1)), -sin(1) + sin(3)/3 - sin(6)/3 + sin(2)), 0.20800494410405)" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intg = sympy.integrate(g(x, y), (y, 1, 3), (x, -2, -1))\n", "\n", "Eq(Eq(Integral(g(x, y), (y, 1, 3), (x, -2, -1)), intg), float(intg), evaluate=False)\n", " # Beachten Sie bitte das letzte Argument.\n", " # Eq in Sympy wird verwendet, um zu prüfen, \n", " # ob zwei Ausdrücke mathematisch gleich sind. \n", " # Hier würde die Auswertung ein 'false ergeben', \n", " # da der analyische Ausdrucke natürlich niemals gleich \n", " # sein wird zu dem numerischen. Das würde uns hier\n", " # aber nicht interessieren.\n", " # Sie können die Auswertung mit evaluate = falsch \n", " # verhindern. \n", " \n", " \n", " " ] }, { "cell_type": "markdown", "id": "baa6c35b", "metadata": {}, "source": [ "Schauen Sie sich bitte einmal die Dokumentation von SymPy an. Dort werden Sie sehen, dass es noch andere Befehle gibt zur Evaluation symbolischer Ausdrücke. Zum Beispiel mit der `N()`-Funktion erhalten Sie ein ähnliches Ergebnis. Hier können Sie auch noch die Anzahl der Stellen angeben, mit denen der Ausdruck numerisch evaluiert werden soll. Wenn Sie kein Argument angeben, erfolgt die Evaluation auf 15 Nachkommastallen, was gerade double-precision entspricht." ] }, { "cell_type": "code", "execution_count": 21, "id": "d976e8b2", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.20800494410404955371428429408811418957133462017301\n", "0.208\n" ] } ], "source": [ "print(sympy.N(intg, 50))\n", "print(sympy.N(intg, 3))" ] }, { "cell_type": "markdown", "id": "23e36a54", "metadata": {}, "source": [ "Abschliessend können wir hier noch das letzte Beispiel aus der Vorlesung am Computer lösen. Hier lautete das zu lösende Integral\n", "$$\n", "\\int_{\\sqrt{\\pi}}^{2\\sqrt{\\pi}}2x\\cos x^2 \\mathrm{d}x\n", "$$\n", "Wir definieren uns das ganz einfach im Folgenden als unsere Funktion und evaluieren dann das bestimmte Integral." ] }, { "cell_type": "code", "execution_count": 22, "id": "9c47b185", "metadata": {}, "outputs": [], "source": [ "def f(x):\n", " return 2 * x * sympy.cos(x**2)" ] }, { "cell_type": "code", "execution_count": 23, "id": "9cdd98c7", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int 2 x \\cos{\\left(x^{2} \\right)}\\, dx = \\sin{\\left(x^{2} \\right)}$" ], "text/plain": [ "Eq(Integral(2*x*cos(x**2), x), sin(x**2))" ] }, "execution_count": 23, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = sympy.integrate(f(x))\n", "\n", "Eq(Integral(f(x), x), intf)" ] }, { "cell_type": "code", "execution_count": 24, "id": "d391716a", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle \\int\\limits_{\\sqrt{\\pi}}^{2 \\sqrt{\\pi}} 2 x \\cos{\\left(x^{2} \\right)}\\, dx = 0$" ], "text/plain": [ "Eq(Integral(2*x*cos(x**2), (x, sqrt(pi), 2*sqrt(pi))), 0)" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "intf = sympy.integrate(f(x), (x, sympy.sqrt(sympy.pi), 2 * sympy.sqrt(sympy.pi)))\n", "\n", "Eq(Integral(f(x), (x, sympy.sqrt(sympy.pi), 2 * sympy.sqrt(sympy.pi))), intf, evaluate=False)" ] }, { "cell_type": "markdown", "id": "f3abd670", "metadata": {}, "source": [ "## Numerisches Lösen von Integralen " ] }, { "cell_type": "markdown", "id": "90ed0747", "metadata": {}, "source": [ "Viele der relevanten Integrale in der Vorlesung werden wir analytisch lösen. Praktisch ist es aber so, dass dies in (ferner) Zukunft nicht immer möglich sein wird. Daher würden wir im Folgenden auch einige Grundzüge der numerischen Integration besprechen. \n", "\n", "Als Erinnerung, das Integral allgemein ist definiert als die Fläche unter einer Funktion.\n", "\n", "Hierfür importieren wir zunächst wieder `NumPy` (unsere go-to Bibliothek für numerisches Rechnen in Python) sowie `matplotlib` zum Visualisieren der Ergebnisse." ] }, { "cell_type": "code", "execution_count": 25, "id": "1b222e20", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "id": "5336112a", "metadata": {}, "source": [ "Nun geben wir das `x`-Achsenintervall an, auf dem wir unsere Funktion auswerten möchten und berechnen die Funktionswerte an den entsprechenden Punkten. Beachten Sie, dass wir hier unsere Funktion nicht explizit definieren, sondern einfach direkt auswerten. Im Prinzip könnte man natürlich auch `def f(x): ...` verwenden, ist hier aber nicht zwingend notwendig." ] }, { "cell_type": "code", "execution_count": 26, "id": "9643d256", "metadata": {}, "outputs": [], "source": [ "x = np.linspace(-2.5, 2.5, 100) # 100 Punkte zwischen -2.5 und 2.5\n", "f = np.exp(-(x**2)) # unsere Funktionswerte" ] }, { "cell_type": "markdown", "id": "3cb9e67b", "metadata": {}, "source": [ "Und schlussendlich visualisieren wir das Ganze." ] }, { "cell_type": "code", "execution_count": 27, "id": "910dfa7d", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAARIAAADQCAYAAADLYO4oAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAkRklEQVR4nO3deXhU9b3H8fc3k0xmJiuEgIBgANl3CItsoriAuGHFtS5d4KFVb6+9t5bW5ba21q0u9bogu2FJiLgAgoJSq6JFCciOLIIFCoQt7JCQ5Hf/mIQbQpZJ5pxZv6/n4Wln5sz3922f/D7POXPO+R0xxqCUUv6ICXYDSqnwp0GilPKbBolSym8aJEopv2mQKKX8pkGilPJbbLAbqKtGjRqZjIyMYLehVNRZuXLlQWNMelWfhV2QZGRkkJeXF+w2lIo6IvKv6j7TQxullN80SJRSfrMtSERkqojsF5H11XwuIvKKiGwTkbUi0suuXpRS9rJzj2Q6MLyGz0cAbcv+jQXesLEXpZSNbPux1RjzuYhk1LDJTUCW8d41uFxEUkWkqTFmr109qeAoKSlh5cqVHDt2DID27dvTokWLIHelrBTMszbNgV0VXu8ue++CIBGRsXj3WmjZsmVAmlP+y8/P57nnnmPGzFkc2J9/3mcDBg5k7Jgx3HvvvYhIkDpUVglmkFT111PlmgbGmInARIDMzExd9yAMzJs3j5/+7OccPnyI+GYdaHTDvTiSG4ExFO7eyMpN/+D+++9n9uzZTJs2jWbNmgW7ZeWHYAbJbqDi/u3FwJ4g9aIsYoxh/PjxPPfccxATS9qI/yCx69XnbeNq0YXk/qM5sfpDPvl0Cp27dGHJ4sX06dMnSF0rfwXz9O984N6yszf9gaP6+0h4M8bw6KOPekPEEUfjW5+4IETKiQhJPa/jovtf4RTxXHX1NaxevTqwDSvL2Hn6Nxv4J9BeRHaLyM9EZJyIjCvbZBGwHdgGTAJ+aVcvKjD+8pe/8PTTT3tDZPQfcLeq/Yx+XMPmpN/+FGeI48phw9i0aVMAOlVWk3BbajEzM9PoJfKh56OPPmLEiBFlIfI/uC/pUafvny3Yw4HZj9CmRVNW5uXh8XjsaVTVm4isNMZkVvWZXtmq/LZnzx7uuvvHEBNL41seq3OIAMQ1aEaDkf/Nd999x0MPPWR9k8pWGiTKLyUlJdx5110UFBSQNuJB3K1717uWO6MHyf1HM3XqVLKzsy3sUtlNg0T5ZfLkyXz+2Wck97mZxC5X+V0vddDduJp35JcPPMDBgwct6FAFggaJqrcDBw7wm988grNJa1KH/sSSmhLjoMG1D3D06DHGjx9vSU1lPw0SVW+PPPIIJ06eJO36/7b06lRnegbudpcxZcoUli9fblldZR8NElUvy5cvZ/r06SRl3oSzkbW3LZzZuZbTP6wh1pPCuF/8gtLSUkvrK+tpkKg6M8bw2/HjiUtsQMrAOy2tfWbnWg68/wyNR/2OlCt/zprVq8nNzbV0DGU9DRJVZ0uXLuXzzz4jsd9oYpxuy+qWh0j6zeNxtexGQqfLiW+cweNPPEFxcbFl4yjraZCoOjHG8Lvf/x5nSjpJPUZYVrdyiACIxJA88G62bd1KVlaWZWMp62mQqDpZsGABeStWkHjZHUhsnCU1qwqRcu62/XE1a8f//OEPFBYWWjKesp4GifKZMYYn//QnnA2bkthlmCU1awoR8N7clzzwbnbv2sXMmTMtGVNZT4NE+eyzzz5jZV4eiX1uQRz+r0BRW4iUc7XqRfxFbXj2uef0DE6I0iBRPnv22WeJS0wlofOVftfyNUSgbMmBPrewdcsWFixY4PfYynoaJMona9eu5aOPPiKh5w3ExMX7VasuIVLO02EQztQmPPPss36NreyhQaJ88te//hWH001iz+v8qlOfEAHvpfMJmTez/J//5Msvv/SrB2U9DRJVq3379pGdk4On61U43En1rlPfECmX2PVqYt1JvPzy3+rdg7KHBomq1aRJkyg+e5akXtfXu4a/IQIQ43Th6Xo1773/Hnv26PK+oUSDRNWouLiYNyZMwN2qF3ENm9erhhUhUi6x53WUlJQwceJEv+ooa2mQqBrNmzePvXv2kNhzZL2+b2WIAMSlXoS7dW/emDCBoqIiv+spa2iQqBq9+uprOFOb4G5T5VKdNbI6RMol9RzJ/vx83nvvPctqKv9okKhqfffdd/zjH5/i6T4ciXHU6bt2hQiAq3Vv4pLS+N9XX7W0rqo/DRJVrSlTpiAxDhK71m0JRTtDBKBw13qKT5/gy2XL2Lp1q+X1Vd3ZGiQiMlxENovINhG5YN08EUkRkQUiskZENoiINev1Kb8VFRUxbfp0XJf2xZHQwOfv2R0i5fXTRj6MxDiYOnWq5WOourPzAVkO4DVgBNAJuFNEOlXa7AFgozGmOzAUeEFEnHb1pHz3wQcfcOjgQZK6XePzdwIVIuk3jyehwyBcrXszddo0XaskBNi5R9IX2GaM2W6MKQJygJsqbWOAJPEu+JkIHAb0ryIETJo0ibjkRrh8eFoeBDZEyusndruW/fn5LFq0yPLxVN3YGSTNgV0VXu8ue6+iV4GOeB8evg74lTHmgts7RWSsiOSJSN6BAwfs6leV2bVrF0uWLMHTeZhPP7IGI0QA3G0yiUtqyOTJUywfU9WNnUFS1bLilZ8Pei2wGmgG9ABeFZHkC75kzERjTKYxJjM9Pd3qPlUlM2bMoLS0lIRuVT8AvKJghQh4779xd7qSRYsWkp+fb/nYynd2BsluoEWF1xfj3fOo6CfAu8ZrG7AD6GBjT6oWxhimv/UWrhZdiEu9qMZtgxki5RK6XElJSYk+mS/I7AySFUBbEWlV9gPqHcD8StvsBIYBiEgToD2w3caeVC1WrFjB1i1b8NSy5kgohAiAs1FLXE3bMv2ttyzvQfnOtiAxxhQDDwKLgU1ArjFmg4iME5FxZZv9CRggIuuApcBvjTH6nMYgysrKIibOSUKHgdVuEyohUs7d+UrWrF7NunXrLO9F+cbW60iMMYuMMe2MMW2MMU+VvTfBGDOh7L/vMcZcY4zpaozpYozRRTmDqKioiFmzZ+O6tD8x8QlVbhNqIQKQ0HEIEuPQleaDSK9sVecsWrSIIwUFJFZzWBOKIQLg8KTgap1J1owZek1JkGiQqHOysmYQl9gAV6ueF3wWqiFSLqHLlezPz2fp0qWW96Zqp0GiADhy5AgLFy7E1X7wBdeOhHqIAHja9MHhStCzN0GiQaIAeO+99ygqKiSh0+XnvR8OIQIgsU5cbS9j7jvvcvr0aQs7VL7QIFEAzJo1C2fDpjibtjv3XriESLmETkM5eeI4CxcutKA7VRcaJIq9e/fy6aef4u5wOd7bnsIvRABcLbsSl9SQmTNnWVJP+U6DRDFnzhzvJfEdvYc14Rgi4L1k3tV+MIsWLaKgoMCyuqp2GiSKmbNmEX9RG+IatQjbECkX17A5Z88W8e6771peW1VPgyTKff/996zMy8PdYUjYh8iZnWsp+HwGccmN9OxNgGmQRLk5c+YAEJucHvYhcuD9Z2g86ne4Ow/j008/1TuCA0iDJMpl5+TgTL+Ewx9PCPsQKa+f0HEwpaWlvPPOO5aPpaqmQRLFNm7cyPp16zh7JD9iQgTAmZ5BfPolzNbDm4DRIIli5Yc1aSP/M2JCpJyr/SC++vJLdu/ebfm46kIaJFHKGMOMmbNwNmlNQvtBltcP9g+3CR2HYIzh7bfftnxsdSENkii1Zs0admz/nsQeIyyvHewQAe9p4PiL2pCdk2P5+OpCGiRRKjc3F4lx4Gk3wNK6oRAi5dztB7Pim2/44YcfLO9DnU+DJAoZY8iZM4f4lt1weFIsqxtKIQLg6eA9ZJs7d67lvajzaZBEoW+//ZYd27fj6TDYspqhFiIAcakX4WrWTg9vAkCDJAr9/2FNf0vqhWKIlHO1G8SqlSvZvl3XFLeTBkmUMcaQPWcOroweONwXPEKozkI5RAASyg5v9OyNvTRIokxeXh47f/gBjwWnfEM9RABiUxrjatZeD29spkESZXJzc4lxxOJud5lfdcIhRMq52w9izerVbNu2zaLuVGW2BomIDBeRzSKyTUTGV7PNUBFZLSIbROQzO/uJdsYY5uTmEp/RA4crsd51wilEADxlz+jRwxv72BYkIuIAXgNGAJ2AO0WkU6VtUoHXgRuNMZ2B0Xb1o7xP0du1cyee9vU/WxNuIQIQm9wYV/OO5JTdEqCsZ+ceSV9gmzFmuzGmCMgBbqq0zV14n/27E8AYs9/GfqJebm4u4ojF07Zfvb4fjiFSzt1+EGvXrGHr1q2W1lVedgZJc2BXhde7y96rqB3QQET+ISIrReReG/uJasYYcnJzcWX0JKYehzXhHCIAjqSGgB7e2MXOIJEq3jOVXscCvYGRwLXA4yLSrvKXRGSsiOSJSN6BAwes7zQKfPPNN/x7165zV3vWRbiHyJmdazm85A2c6Zfo4Y1N7AyS3UCLCq8vBvZUsc1HxpiTZQ8P/xzoXrmQMWaiMSbTGJOZnp5uW8ORLDc3l5jYODyX1u2wJhJCpLx+QrdrWLd2LZs3b7Z8nGhnZ5CsANqKSCsRcQJ3APMrbTMPGCwisSLiAfoBm2zsKSqVlpaSk5tL/CU96nRYE0kh4mrZDU97PXtjF9uCxBhTDDwILMYbDrnGmA0iMk5ExpVtswn4CFgLfANMNsast6unaPX111+zZ/duPB2H+PydSAsRgNikRrgu7ky2Ht5YLtbO4saYRcCiSu9NqPT6eeB5O/uIdt7DGqfPhzWRGCLl3B0GsfGTN9m0aRMdO3a0fOxo5dMeiYg0FpFRIvKAiPxURPqKiF4VGwbOHda06kVMvKfW7SM5RADv+isi5ObmWj52NKsxDETkChFZDCzEe2FZU7wXlz0GrBORP4qI/3d+Kdt89dVX7Nuz59zNazUJ9iQPRP3YpDRcF3fWszcWq+3Q5jpgTPkFYxWJSCxwPXA1oOv+h6jc3Fxi4py42/StcbtQmOSBqu/pMIjvPp7Ahg0b6Ny5s+W9RKMa90iMMb+pKkTKPis2xrxvjNEQCVElJSXMeftt4lv1rvGwJpQmeSDqe9oNRCTm3Cr6yn++/kYyQ0RSKrzOEJGl9rWlrPDFF1+wf98+EmpYCS3UJnkg6jsSGxDfsgvZOTkYU/kaSVUfvv5gugz4WkSuE5ExwBLgZdu6UpbIycnB4XRVe1gTipM8UPU9HYawbetW1qxZY3lf0cinIDHGvAn8HO8FZE8CQ4wxC+xsTPnn7NmzvD13LvFt+hLjdF3weShP8kDU97S7DIlxkKMLHlnC10Obe4CpwL3AdGCRiFxwKbsKHX//+985fOhQlYc1oT7JA1Hf4UnBldGD2Xp4YwlfD21+BAwyxmQbY34HjMMbKCpEzZkzh1hXAu7Wvc97PxwmeaDqezoMYde//sWKFSss6i56+Xpoc3PFtUKMMd/gvS9GhaDCwkLeefc94i/th8Q6z70fTpM8EPU97foTExunhzcWqO2CtMdEpGFVnxljikTkShG53p7WVH0tXryYY0ePkNDh/++tCbdJHoj6MfEJxLfqTXZODiUlJZbUjFa1XZC2DlggImeAVcABwAW0BXoAnwB/sbNBVXfZ2dnEepJxZfQAwnOSB6p+Qsch7Jv/HF988QVDhw61tHY0qe3Q5lZjzEC8d/BuABzAMWAm0NcY87AxRlcaCiEnTpzg/ffn4Wo/CHHEhvUkD0T9mHgPSAzZ2dmW144mtQVJbxG5BLgb71oibwJZeNcacdvcm6qH+fPnc+bMaRI6XR72kzwQ9Q9+8CKult2Yk/s2RUVFlo8RLWoLkgl41wvpAORV+Ley7D9ViJk1azbOlHRMydmwn+SBqp/U5yaOHilgyZIllo8TLWq71+YVY0xHYKoxpnWFf62MMa0D1KPy0aFDh1i8ZDGxzTpycN5zETHJA1HfndGTWE8ys2fPtnysaOHr6d9f2N2I8t/cuXMpKS7m9PaVETPJA1FfHLG42g3gvffnceLECcvHjAa6OFEEycqaATGxpI/6fcRM8kDVT+h8BWdOn+L999+3fNxooEESIXbs2MFXX31JYtdhuC+x/u6FSA4RgPjmHXGmNiErK8vysaOBBkmEmDlzJgApA263vHawJ3kg6ovE4Oo4lKVLl7J3717Le4h0GiQRwBjDW1kzcLXsSmxyY0trh8IkD1T9xM5DvWvc6iXzdaZBEgFWrFjB99u24ul0haV1Q2mSB6J+XFoLXE3b8pYe3tSZrUEiIsNFZLOIbBOR8TVs10dESkTkVjv7iVQzZswgJtZJQoeBltUMtUkeqPruTlewZvVqNmzYYHlPkcy2IBERB/Aa3tXnOwF3ikinarZ7Fu9l+KqOCgsLmTlrFq5L+xETn2BJzVCd5IGon9BxCBLj4K233rK8r0hm5x5JX2CbMWa7MaYIyAFuqmK7h/CuQr+/is9ULRYuXMiRggISu15lSb1QnuSBqO9ISMXVpg9vZWVRXFxseX+Rys4gaQ7sqvB6d9l754hIc2AU3kvxVT1MmzaduOS0c3f6+iPUJ3mg6id2Gcb+/HwWL9adZF/ZGSRSxXuV17R7GfitMabGxSBEZKyI5IlI3oEDerNxufz8fD78cBHujlcgMQ6/aoXLJA9EfXebTGITUpg+fbo1zUUBO4NkN9CiwuuLgT2VtskEckTkB+BW4HURublyIWPMRGNMpjEmMz093aZ2w8/MmTMpKSkhseswv+qE0yQPRH1xxOHucDnz5s/n0KFDFnQY+ewMkhVAWxFpJSJO4A68SxGcU3bzX4YxJgOYC/zSGPO+jT1FDGMMU6dNw9W8A3FpLWr/QjXCbZIHqn5it6s4W1SkN/L5yLYgMcYUAw/iPRuzCcg1xmwQkXEiMs6ucaPF119/zcYNG/B0qf+PrOE6yQNR39m4Nc6GzZnw5pu6yrwPaltq0S/GmEXAokrvVfnDqjHmfjt7iTSTJk3C4XST0HFI7RtXIZwneaDqnz1+iI2H/82KFSvo27fmZydHO72yNQwdO3aM2dk5uDsMrvGZvtWJhEkeiPqNbnwEh9PFpEmTLB8j0miQhKHs7GzOnD5FYvdr6/zdSJnkgajvubQvrvaDmDU7m+PHj1s+ViTRIAlDb06cSHzjDJxN29Xpe5E0yQNVP6n7tZw+dVJv5KuFBkmYWbVqFd+uWoWn27WIVHWpTtUicZIHor6zWQfi0y/hjQl6zWRNNEjCzKuvvorD6SKxs+93+kbqJA9EfRHB0304365axTfffGP52JFCgySMHD58mNmzs3F3vJwYV6JP34nkSR6o+oldhuGId/P6669bPn6k0CAJI9OmTaOw8AxJvUb6tH0oTMJIqB8T78Hd6Qqys3M4ePCg5X1EAg2SMFFaWsprr7+O6+JOOBvX/iSQUJmEkVI/qed1FBUVMnXqVMt7iQQaJGFiyZIl7Ni+ncSete+NhNokjIT6zvQMXC268Nrrr+sDx6ugQRImXnrpJeKSGuJpP6DG7UJxEkZK/cRe17PzX/9iwYIFlvcV7jRIwsD69etZsmQJCT1GIo64arcL5UkYCfU97S7DmdqEv77wguW9hTsNkjDw0ksv4YiLJ7HH8Gq3CfVJGAn1JcZBQq8b+HLZMvLy9NHXFWmQhLj9+/czc+Ys3J2vwOFJqXKbcJiEkVI/sds1OOI9vPTSSxZ1Fxk0SELcG2+8QVFRIcmZVS13G16TMBLqx8R78HS7hjm5uezatav2L0QJDZIQduLECV7+2yu42/SpcvGicJuEkVI/ufeNlBrDiy++aFnNcKdBEsImTZrEkYLDpFx22wWfheskjIT6sSmNiW/ZnTfemICuIeylQRKiCgsLee7553G37Ep8847nfRbOkzBS6hfu2UxhUSGvvPKK5fXDkQZJiMrKymLf3r0k9R993vuRMAkjoX7jWx7F0/Yy/vbK/3Ls2DHLxwk3GiQh6OzZszz9zDPEN70UV0bPc+9HyiSMlPrJ/Udz/NhRXnvtNcvHCjcaJCEoKyuLHdu3k3zZHefWHIm0SRgJ9eObtsXdJpNnn3ueo0ePWj5mONEgCTGFhYX84Y9/JL5pO9yX9gMicxJGSv2UQT/m6JECXn75ZcvHDScaJCFm8uTJ7N61i5TBP0ZEInoSRkL9+IsuxdNuAH994UUOHz5s+fjhwtYgEZHhIrJZRLaJyPgqPr9bRNaW/ftKRLrb2U+oO336NH/6859xteiMK6Nn0CeJ1vetfsqguzhx4jjPP/+85T2EC9uCREQcwGvACKATcKeIdKq02Q7gcmNMN+BPwES7+gkHL730Evn79pEy+B4Kd60LiUmi9Wuv70zPwNNxCC+//LeovdrVzj2SvsA2Y8x2Y0wRkAOcd523MeYrY0xB2cvleJ8PHJX27dvHU395Gk/b/mBKQ2aSaH3f6jcYch9FxSU8+uijlvcTDuwMkuZAxXjeXfZedX4GfGhjPyHtiSee4PTp07jbDQi5SaL1axeb0pjE3jcyY8aMqLwz2M4gqepZCVU+RFVErsAbJL+t5vOxIpInInmReEny+vXrmTJlCq5L+1Lw98khN0m0vm9SLhtNXEIqD//611H3vGA7g2Q3UPFOs4uBPZU3EpFuwGTgJmPMoaoKGWMmGmMyjTGZ6enptjQbLMYYfvnAA0icizM714XsJNH6tYuJTyBp0N0s++KLqHuglp1BsgJoKyKtRMQJ3AHMr7iBiLQE3gXuMcZssbGXkDVjxgy++PxzSktLaTzq9yE7SbS+bxK7XYOraVt+9fDDUXWRmm1BYowpBh4EFgObgFxjzAYRGSci48o2ewJIA14XkdUiElUHlwUFBTz86/+CmFjSf/RYyE8SrV87iXGQes0DHNx/gMcee8yCDsODhNuxXGZmpomUH7PGjh3LpEmTaDj8QZK6V7+MYn2F2ySMpPqHP57AidWL+Hr5cvr06WNp7WARkZXGmMyqPtMrW4NkyZIlTJo0CU+HwRoiEVjf1aoXIPz4nns5c+aM5fVDjQZJEBw5coR777ufuJTGNBr5sOX1w30SRkL9Q4teJmXIPWzZ/B2PP/645WOEGg2SIPjVr37F/vx80m4aj8Q6La0dCZMwUuqn9LuVxO7DeeGFF/jyyy8tHyuUaJAE2KxZs8jKyiK5/2jim7aztHYkTcJIqd/gip8Sl9qE2++8M6Jv6tMgCaDNmzczZuxYXBd3JmXQXZbWjsRJGAn1Y+I9NLzhEfbu2ct9990XsReqaZAEyKlTp7jlR7dSLHGk3fAbJMZhWe1InYSRUj++aTtShv6UDz74IGJXntcgCYDS0lLuu+8+Nm7cQIPrfk1sciPLagd7kmh93+on9b4BT7sBPPLIIyxevNjyPoJNgyQA/vjHPzJ37lxSL/8J7ta9LasbKpNE69deX0RIG/kwcekZ3Dr6NjZu3Gh5P8GkQWKzGTNm8OSTT5LQ9SqS+46yrG4oTRKt71v9GKebRrc8RqFxcN3Ikezbt8/yvoJFg8RG8+bN4/6f/AT3Jd1Ju/aBcws5+ysUJ4nW901scmPSRj3G7j35DLvqqog5k6NBYpNPPvmE0bfdhrPJpTQa9SjiiLOkbihPEq3vm/hm7Ukb9Sjfbd7C8BEjOH78uMVdBp4GiQ0++OADRo68npjUZjS69Q/ExHssqRsOk0Tr+8ad0YO0G39LXt5KrrjySg4dqnIFjbChQWKx7Oxsbh41ipi0lqTf/hQOd5IldcNpkmh933ja9qfRqN/z7eq1DB5yOXv37rWkbjBokFjEGMOf//xn7rrrLuKadvCGiCfFktrhOEm0vm88l/YjZfA9bNq0id6Zffj2228trR8oGiQWOHnyJHfeeSePP/44CZ2voPFtT+rhjNb3uf7R5W/T8NpfcuhkEQMGDCQ3N9fyceymQeKnVatW0b1HT+bk5pI69H7SRv7ashvxImGSaH3f6id1H07je16ERhncfvvtjBkzhpMnT1o+pl00SOqpqKiIp556in79+rNzfwFN7niKlH636ilerV/v+o6EBjS+4y8k97uVyVOm0LNXL5YtW2b52HbQIKmHzz77jC5du/HYY4/hbNOPJve/YukfWyROEq3vW31xxNFg6P00vu1P/Cv/CIMHD2bMmDGE+tMTNEjqYP369Yy8/nqGDh3Kjj37Sb/1f0i/eTwOd7JlY0TyJNH6vtd3Z/SgyU9fJ7nvLUyZOo1Wrdvw1FNPhezhjq7ZWgtjDMuXL+eZZ55l/vx5IDEkdruaBsPGEBPnsnSsUPkj1vqhVf/swV0UfP4Wp7cup0HDhvzHQw/x4IMP0qiRdTd/+qKmNVs1SKpRUFBATk4OE958k7Vr1uBwJVBSfJb0Gx/xPlbTYqH6R6z1Q6d+4b+/4+jytzm97WvinE5u/dGPGDNmDEOGDMHhsG5ZiupokPjo3//+Nx999BHvvPMuH3/yMcVnzxLfpDXOFl04sf5TGo/6Xcj+kWn96KlfdHAnJ75dxOmNn1J85iSNmzTh9ttu44YbbmDw4MG4XNbuKZcLWpCIyHDgb4ADmGyMeabS51L2+XXAKeB+Y8yqmmpaFSTGGL7//ntWrFjBsmXL+PQf/2BT2a3dzpTGuNoPwtNhMKVFpzg479mw+SPT+tFTv/TsGY4tf5ujy98hJkYoLT6Ly+Vm0KBBDBkymIEDB9KrVy9SU1MtGa+mIIm1ZISqB3UArwFX43185woRmW+MqbgQwwigbdm/fsAbZf9pCWMMR44cYdeuXezYsYPvv/+ezZs3s2HDRtauW8fxY94noTni3cQ160Dq0Ptxt84krtEliAhndq7VENH6IVu/aO8Wjn/7IU1ufxLnRe0o3LWO09tX8sXa9Xyy9BMo20loeUkG3bt1pWPHjrRt25bWrVuTkZFB8+bNiY+Pt6QX24IE6AtsM8ZsBxCRHOAmoGKQ3ARkGe9u0XIRSRWRpsaYam86OH78OAsWLODkyZOcPHmS48ePc/ToUY4ePUpBQQGHDh3iwIGD7N23j/z8fM6cPnXe92M9ycSmtSC29QAaNmlDfNO2xKVnXLD0Ybj/kWn96KvvbtMHdxvvw7hKzpygaO8WivZt4+D+HSz5eh0fLPoQU1J8Xp0GDdNo0qQJzZpeRFpaGmlpaaSkpJCamkpSUhKJiYkkJibi8dR8pbadQdIc2FXh9W4u3NuoapvmQLVBsmXLFm688cYL3o9zuXF6konzJOPwJBGf1oamGX1wpabjSk3Hk9YUT6PmxHmSar1o7OjOzWxaOZcBD7xASsv2Nf+vrAetr/Xtr58MbZoBQ8+9U1pSzJkjBzh1cA+nC/I5c+QAZ44eoOB4Afk/HKBk4/ecPXmUwlMnKK0UOLWx7TcSERkNXGuM+XnZ63uAvsaYhypssxB42hizrOz1UuARY8zKSrXGAmPLXrYHNtvSNDQCDtpU2w7h1i+EX8/h1i/Y1/Mlxpj0qj6wc49kN9CiwuuLgT312AZjzERgotUNViYiedX9mBSKwq1fCL+ew61fCE7Pdl7ZugJoKyKtRMQJ3AHMr7TNfOBe8eoPHK3p9xGlVGiybY/EGFMsIg8Ci/Ge/p1qjNkgIuPKPp8ALMJ76ncb3tO/P7GrH6WUfew8tMEYswhvWFR8b0KF/26AB+zsoY5sP3yyWLj1C+HXc7j1C0HoOeyubFVKhR69+1cp5TcNkgpE5HkR+U5E1orIeyKSGuyeaiMio0Vkg4iUikjInl0QkeEisllEtonI+GD3UxsRmSoi+0VkfbB78ZWItBCRT0VkU9nfxK8CNbYGyfk+BroYY7oBW4DfBbkfX6wHbgE+D3Yj1alwu8QIoBNwp4h0Cm5XtZoODA92E3VUDPyXMaYj0B94IFD/P2uQVGCMWWKMKb+kbzne61pCmjFmkzHGrgv0rHLudgljTBFQfrtEyDLGfA6E1WPwjDF7y296NcYcBzbhvVLcdhok1fsp8GGwm4gQ1d0KoWwiIhlAT+DrQIxn6+nfUCQinwAXVfHRo8aYeWXbPIp3N3FWIHurji89h7iqbm7S04U2EZFE4B3gP40xxwIxZtQFiTHmqpo+F5H7gOuBYSZEzo3X1nMY8OlWCOU/EYnDGyKzjDHvBmpcPbSpoGwhpt8CNxpjTtW2vfKZL7dLKD+VLRQ2BdhkjHkxkGNrkJzvVSAJ+FhEVovIhNq+EGwiMkpEdgOXAQtFZHGwe6qs7Afs8tslNgG5xpgNwe2qZiKSDfwTaC8iu0XkZ8HuyQcDgXuAK8v+fleLyHWBGFivbFVK+U33SJRSftMgUUr5TYNEKeU3DRKllN80SJRSftMgUUr5TYNEKeU3DRJlOxHpU7bGi0tEEsrWyugS7L6UdfSCNBUQIvJnwAW4gd3GmKeD3JKykAaJCoiye2xWAGeAAcaYkiC3pCykhzYqUBoCiXjvZXIFuRdlMd0jUQEhIvPxrozWCmhqjHkwyC0pC0XdeiQq8ETkXqDYGDO7bP3Wr0TkSmPM34Pdm7KG7pEopfymv5EopfymQaKU8psGiVLKbxokSim/aZAopfymQaKU8psGiVLKbxokSim//R9fCbyuLT5CRQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1, figsize=(4, 3)) # wir initialisieren die Grafik\n", "ax.plot(x, f, color=\"black\") # x und f sind gleich große Arrays, die Punkt für Punkt geplottet werden\n", "ax.fill_between(x, f, hatch=\"/\") # die Fläche unter der Funktion\n", "ax.set_xlim(x[0], x[-1]) # x-Achsen limits\n", "ax.set_ylim(0, None) # y-Achsen limits\n", "ax.set_xlabel(\"x\") # x-Achsenbeschriftung\n", "ax.set_ylabel(\"f(x)\") # x-Achsenbeschriftung\n", "plt.show() # und nun zeigen wir das Ganze an" ] }, { "cell_type": "markdown", "id": "d1395274", "metadata": {}, "source": [ "Numerisch werden wir keine unbestimmten Integrale lösen können sondern nur bestimmte. Daher betrachten wir im Folgenden numerische Integrationsverfahren zur Berechnung des bestimmten Integrals einer Funktion $f(x)$ auf einem Intervall $[a,b]$. \n", "$$\n", "I=\\int_a^b f(x)\\mathrm{d}x\n", "$$\n", "Hierfür definieren wir uns zunächst ein beliebiges Integrationsinterval." ] }, { "cell_type": "code", "execution_count": 28, "id": "8654e823", "metadata": {}, "outputs": [], "source": [ "a, b = -1.2, 0.7" ] }, { "cell_type": "markdown", "id": "717406f2", "metadata": {}, "source": [ "Und dann sieht das Ganze folgerndermaßen aus." ] }, { "cell_type": "code", "execution_count": 29, "id": "f1ce6795", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAARIAAADQCAYAAADLYO4oAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAlNklEQVR4nO3deXQUdbr/8feTjayELSB7AFkCQVkCIosy4AbIGBGRRRgcheEqOt7Ri3BFHbeRxVHHEWEA2QKEBEc2iYAig+CaoOwIsv0gsmPYwhIIz++PTrghBLJ0VXdX8n2dw9Hurn7q0cP3c6q6qr5fUVUMwzDc4eftBgzDcD4TJIZhuM0EiWEYbjNBYhiG20yQGIbhNhMkhmG4LcDbDRRXlSpVNDo62tttGEaZs27dumOqGlXQZ44LkujoaNLS0rzdhmGUOSLy/673mTm1MQzDbSZIDMNwm21BIiLTROSIiGy+zuciIu+LyE4R2SgirezqxTAMe9l5RDIDuO8Gn3cDGub8GQpMtLEXwzBsZNuPrar6lYhE32CTB4BZ6npq8DsRqSAi1VX1oF09Gd6RnZ3NunXrOHXqFACNGzemdu3aXu7KsJI3r9rUBPbneZ2e8941QSIiQ3EdtVCnTh2PNGe47/Dhw4wbN46E2XM4euTwVZ+179CBoUOGMGjQIETESx0aVvFmkBT0t6fAOQ1UdTIwGSAuLs7Me+AAixYt4o+PP8Fvvx2nXI0mVOk5CP/yVUCVC+lbWbftPwwePJi5c+cyffp0atSo4e2WDTd4M0jSgbzHt7WAA17qxbCIqjJy5EjGjRsHfgFU7vYM4c3vvmqb4NqxlG/3MGfWf8YXqz6iWWwsK5Yvp02bNl7q2nCXNy//LgYG5Vy9aQecNL+POJuq8uKLL7pCxD+Qqr1fviZEcokIES27c9Pg9zlLOe66+x7Wr1/v2YYNy9h5+TcR+BZoLCLpIvK4iAwTkWE5m6QAu4GdwBTgSbt6MTzjb3/7G2+99ZYrRB7+KyH1Cr+iH1ipJlGPvMl5AunStSvbtm3zQKeG1cRpUy3GxcWpuUXe9yxbtoxu3brlhMgrhNRtUazvX8w4wNG5I2hQuzrr0tIIDQ21p1GjxERknarGFfSZubPVcNuBAwfoP+BR8Augaq/RxQ4RgMCKNajY43l+/vlnnn76aeubNGxlgsRwS3Z2Nv369ycjI4PK3YYTUr91iWuFRLegfLuHmTZtGomJiRZ2adjNBInhlqlTp/LV6tWUbxNPeOxdbter0HEAwTVjePKppzh27JgFHRqeYILEKLGjR4/yP/8zgqBq9anQ+TFLaoqfPxXvfYqTJ08xcuRIS2oa9jNBYpTYiBEjOJOZSeX7n7f07tSgqGhCGt3ORx99xHfffWdZXcM+JkiMEvnuu++YMWMGEXEPEFTF2scWzu/byLm9GwgIjWTYf/0Xly9ftrS+YT0TJEaxqSovjBxJYHhFIjv0s7T2+X0bObpwDFUfHEVklyfYsH49ycnJlu7DsJ4JEqPYVq5cyVerVxN+28P4BYVYVjc3RKLiRxJc5xbCmt5JuarRvPTyy1y6dMmy/RjWM0FiFIuqMup//5egyCgiWnSzrG7+EAEQ8aN8hwHs/OUXZs2aZdm+DOuZIDGKZcmSJaSlphJ+e18kINCSmgWFSK6Qhu0IrtGIV/76Vy5cuGDJ/gzrmSAxikxVee311wmqVJ3w2K6W1LxRiIDr4b7yHQaQvn8/s2fPtmSfhvVMkBhFtnr1atalpRHephfi7/4MFIWFSK7geq0od1MDxo4bZ67g+CgTJEaRjR07lsDwCoQ16+J2raKGCORMOdCmF7/s2MGSJUvc3rdhPRMkRpFs3LiRZcuWEdayJ36B5dyqVZwQyRXapCNBFaoxZuxYt/Zt2MMEiVEkb7/9Nv5BIYS37O5WnZKECLhunQ+Li+e7b7/l66+/dqsHw3omSIxCHTp0iMR58whtfhf+IRElrlPSEMkV3vxuAkIieO+9f5S4B8MeJkiMQk2ZMoVLFy8S0er+EtdwN0QA/IKCCW1+NwsWLuDAATO9ry8xQWLc0KVLl5g4aRIh9VoRWKlmiWpYESK5wlt2Jzs7m8mTJ7tVx7CWCRLjhhYtWsTBAwcIb9mjRN+3MkQAAivcREj91kycNImsrCy36xnWMEFi3NAHH0wgqEI1QhoUOFXnDVkdIrkiWvbgyOHDLFiwwLKahntMkBjX9fPPP/Of/6wi9Nb7ED//Yn3XrhABCK7fmsCIyvzzgw8srWuUnAkS47o++ugjxM+f8ObFm0LRzhABuLB/M5fOneHrtWv55ZdfLK9vFJ+tQSIi94nIdhHZKSLXzJsnIpEiskRENojIFhGxZr4+w21ZWVlMnzGD4Jvb4h9WscjfsztEcutX7vHfiJ8/06ZNs3wfRvHZuUCWPzAB6AY0BfqJSNN8mz0FbFXVW4HOwN9FJMiunoyi+/TTTzl+7BgRt9xT5O94KkSi4kcS1qQjwfVbM236dDNXiQ+w84ikLbBTVXerahYwD3gg3zYKRIhrws9w4DfA/K3wAVOmTCGwfBWCi7BaHng2RHLrh99yL0cOHyYlJcXy/RnFY2eQ1AT253mdnvNeXh8AMbgWD98E/FlVr3m8U0SGikiaiKQdPXrUrn6NHPv372fFihWENutapB9ZvREiACEN4giMqMTUqR9Zvk+jeOwMkoKmFc+/Pui9wHqgBtAC+EBEyl/zJdXJqhqnqnFRUVFW92nkk5CQwOXLlwm7peAFwPPyVoiA6/mbkKZdSElZyuHDhy3ft1F0dgZJOlA7z+tauI488noM+ERddgJ7gCY29mQUQlWZMXMmwbVjCaxw0w239WaI5AqL7UJ2drZZmc/L7AySVKChiNTL+QG1L7A43zb7gK4AIlINaAzstrEnoxCpqan8smMHoYXMOeILIQIQVKUOwdUbMmPmTMt7MIrOtiBR1UvAcGA5sA1IVtUtIjJMRIblbPY60F5ENgErgRdU1azT6EWzZs3CLzCIsCYdrruNr4RIrpBmXdiwfj2bNm2yvBejaGy9j0RVU1S1kao2UNU3c96bpKqTcv79gKreo6rNVTVWVc2knF6UlZXFnLlzCb65HX7lwgrcxtdCBCAs5g7Ez9/MNO9F5s5W44qUlBROZGQQfp3TGl8MEQD/0EiC68cxKyHB3FPiJSZIjCtmzUogMLwiwfVaXvOZr4ZIrrDYLhw5fJiVK1da3ptROBMkBgAnTpxg6dKlBDfudM29I74eIgChDdrgHxxmrt54iQkSA4AFCxaQlXWBsKZ3XvW+E0IEQAKCCG54Ox//+xPOnTtnYYdGUZggMQCYM2cOQZWqE1S90ZX3nBIiucKadibzzGmWLl1qQXdGcZggMTh48CCrVq0ipMmduB57cl6IAATXaU5gRCVmz55jST2j6EyQGCQlJbluiY9xndY4MUTAdct8cONOpKSkkJGRYVldo3AmSAxmz5lDuZsaEFiltmNDJFdgpZpcvJjFJ598Ynlt4/pMkJRxu3btYl1aGiFN7nB8iJzft5GMrxIILF/FXL3xMBMkZVxSUhIAAeWjHB8iRxeOoeqDowhp1pVVq1aZJ4I9yARJGZc4bx5BUXX57fNJjg+R3PphMZ24fPky//73vy3fl1EwEyRl2NatW9m8aRMXTxwuNSECEBQVTbmousw1pzceY4KkDMs9ranc49lSEyK5ght35JuvvyY9Pd3y/RrXMkFSRqkqCbPnEFStPmGNO1pe39s/3IbF3IGqMn/+fMv3bVzLBEkZtWHDBvbs3kV4i26W1/Z2iIDrMnC5mxqQOG+e5fs3rmWCpIxKTk5G/PwJbdTe0rq+ECK5Qhp3IvWHH9i7d6/lfRhXM0FSBqkq85KSKFfnFvxDIy2r60shAhDaxHXK9vHHH1vei3E1EyRl0E8//cSe3bsJbdLJspq+FiIAgRVuIrhGI3N64wEmSMqg/zutaWdJPV8MkVzBjTry47p17N5t5hS3kwmSMkZVSUxKIji6Bf4h1ywhVGy+HCIAYTmnN+bqjb1MkJQxaWlp7Nu7l1ALLvn6eogABERWJbhGY3N6YzMTJGVMcnIyfv4BhDS63a06TgiRXCGNO7Jh/Xp27txpUXdGfrYGiYjcJyLbRWSniIy8zjadRWS9iGwRkdV29lPWqSpJycmUi26Bf3B4ies4KUQAQnPW6DGnN/axLUhExB+YAHQDmgL9RKRpvm0qAB8Cv1fVZsDDdvVjuFbR279vH6GNS361xmkhAhBQvirNW7UhOTnZknrGtew8ImkL7FTV3aqaBcwDHsi3TX9ca//uA1DVIzb2U+YlJycj/gGENrytRN93Yojk6to9nvXr1/PLL79YWtdwsTNIagL787xOz3kvr0ZARRH5j4isE5FBNvZTpqkq85KTCY5uiV8JTmucHCIA7Tr9DjCnN3axM0ikgPc03+sAoDXQA7gXeElEGuX/kogMFZE0EUk7evSo9Z2WAT/88AO/7t9/5W7P4nB6iJzft5EzZ07Rvn17c3pjEzuDJB2oned1LeBAAdssU9XMnMXDvwJuzV9IVSerapyqxkVFRdnWcGmWnJyMX0AgoTcX77SmNITI0YVjqNegMX369GHDhg1s377d8v2UdXYGSSrQUETqiUgQ0BdYnG+bRUAnEQkQkVDgNmCbjT2VSZcvX2ZecjLl6rYo1mlNaQmRqPiRhJcvT+/evQFzemMH24JEVS8Bw4HluMIhWVW3iMgwERmWs802YBmwEfgBmKqqm+3qqaz6/vvvOZCeTmjMHUX+TmkKkdz6NWvWpGPHjub0xgYBdhZX1RQgJd97k/K9Hg+Mt7OPss51WhNU5NOa0hgiufr06cMzzzzDtm3biImJsXzfZVWRjkhEpKqIPCgiT4nIH0WkrYiYu2Id4MppTb1W+JULLXT70hwiAA899BAiYo5KLHbDMBCR34nIcmAprhvLquO6uWw0sElEXhUR95/8MmzzzTffcOjAgSsPr92Itwe5J+rXqFGDTp06mSCxWGGnNt2BIbk3jOUlIgHA/cDdgJn330clJyfjFxhESIO2N9zOFwa5p+r36dOH4cOHs2XLFpo1a2Z5L2XRDY9IVPV/CgqRnM8uqepCVTUh4qOys7NJmj+fcvVa3/C0xpcGuSfqP/TQQ/j5+V2ZRd9wX1F/I0kQkcg8r6NFZKV9bRlWWLNmDUcOHSLsBjOh+dog90T9m266iTvvvJOkpCRU898jaZREUX8wXQt8LyLdRWQIsAJ4z7auDEvMmzcP/6Dg657W+OIg91T9vn37smPHDjZs2GB5X2VRkYJEVf8FPIHrBrLXgDtUdYmdjRnuuXjxIvM//phyDdriFxR8zee+PMg9Ub9Xr14EBAQwz0x4ZImintoMBKYBg4AZQIqIXHMru+E7vvzyS347frzA0xpfH+SeqF+lShXuuusuc3pjkaKe2jwEdFTVRFUdBQzDFSiGj0pKSiIgOIyQ+q2vet8Jg9xT9fv27cvevXtJTU21qLuyq6inNvF55wpR1R9wPRdj+KALFy7w708WUO7m25CAoCvvO2mQe6J+fHw8QUFB5vTGAoXdkDZaRCoV9JmqZolIFxG5357WjJJavnw5p06eIKzJ/z1b47RB7on6kZGRdOvWjaSkJLKzsy2pWVYVdkSyCVgiIitFZLyIjBCRl3MuB28CegLf29+mURyJiYkEhJYnOLoF4MxB7qn6/fr148CBA6xZs8bSumVNYUHSW1U74HqCdwvgD5wCZgNtVfW/VdXMNORDzpw5w8KFiwhu3BHxD3D0IPdE/Z49exIWFkZiYqLltcuSwoKktYjUBQbgmkvkX8AsXHONhNjcm1ECixcv5vz5c4Q1vdPxg9zu+gChoaHEx8czf/58srKybNlHWVBYkEzCNV9IEyAtz591Of80fMycOXMJioxCsy86epB7IkRy9evXj4yMDFasWGHrfkqzwp61eV9VY4Bpqlo/z596qlrfQz0aRXT8+HGWr1hOQI0Yji0a59hB7skQAbjnnnuoXLkyc+fOtX1fpVVRL//+l92NGO77+OOPyb50iXO71zl2kHs6RAACAwPp3bs3ixYt4syZMx7ZZ2ljJicqRWbNSgC/AKIe/F9HDnJvhEiuRx99lLNnz7Jw4UKP7re0MEFSSuzZs4dvvvma8OZdCalr/dMLTg+RM6dO3fDz9u3bEx0dTUJCguX7LgtMkJQSs2fPBiCy/SOW13Z6iJzft5E9u268BIWfnx8DBgzgiy++4ODBg5b3UNqZICkFVJWZsxIIrtOcgPJVLa1dGkIkd12bwjz66KOuOW7NLfPFZoKkFEhNTWXXzl8Ibfo7S+uWlhDJXdemME2aNCEuLs6c3pSArUEiIveJyHYR2SkiI2+wXRsRyRaR3nb2U1olJCTgFxBEWJMOltUsTSFSnPoDBw7kp59+YsuWLZb3VJrZFiQi4g9MwDX7fFOgn4g0vc52Y3Hdhm8U04ULF5g9Zw7BN9+GX7kwS2r66iD3RP2+ffsSEBDAzJkzLe+rNLPziKQtsFNVd6tqFjAPeKCA7Z7GNQv9kQI+MwqxdOlSTmRkEN78Lkvq+fIg90T9qlWr0qNHDxISErh06ZLl/ZVWdgZJTWB/ntfpOe9dISI1gQdx3YpvlMD06TMILF/5ypO+7vD1Qe6p+oMHD+bQoUMsX24OkovKziCRAt7LP6fde8ALqnrDySBEZKiIpIlI2tGj5mHjXIcPH+azz1IIifkd4ufvVi2nDHJP1O/evTtVqlRhxowZ1jRXBtgZJOlA7TyvawEH8m0TB8wTkb1Ab+BDEYnPX0hVJ6tqnKrGRUVF2dSu88yePZvs7GzCm3d1q46TBrkn6gcFBTFgwAAWL17M8ePHLeiw9LMzSFKBhiJST0SCgL64piK4Iufhv2hVjQY+Bp5U1YU29lRqqCrTpk8nuGYTAivXLvwL1+G0Qe6p+o899hhZWVnmQb4isi1IVPUSMBzX1ZhtQLKqbhGRYSIyzK79lhXff/89W7dsITS25D+yOnWQe6L+rbfeSsuWLZk6daqZZb4IbL2PRFVTVLWRqjZQ1Tdz3pukqtf8uKqqg1X1Yzv7KU2mTJmCf1AIYTF3FL5xAZw8yD1RH2DIkCFs3LjRzDJfBObOVgc6deoUcxPnEdKk0w3X9L0epw9yTz0l3L9/f0JDQ5kyZYpt+ygtTJA4UGJiIufPnSX81nuL/V2nD3JPTjUQGRnJI488QmJiIqdPn7Z1X05ngsSB/jV5MuWqRhNUvVGxvuf0Qe6N+UqGDBlCZmameZCvECZIHObHH3/kpx9/JPSWexEp6Fadgjl9kHtr0qN27drRrFkz/vWvf3lsn05kgsRhPvjgA/yDgglvVvQnfZ0+yL05c5qIMGzYMNatW8cPP/zg0X07iQkSB/ntt9+YOzeRkJg78QsOL9J3nD7IvRkiuQYNGkR4eDgffvihV/bvBCZIHGT69OlcuHCeiFY9irS90we5L4QIQPny5Rk4cCDz5s3j2LFjXuvDl5kgcYjLly8z4cMPCa7VlKCqha8E4vRB7ishkuvJJ5/kwoULTJs2zdut+CQTJA6xYsUK9uzeTXjLwo9GnD7IfS1EAGJjY7njjjuYOHGiWXC8ACZIHOLdd98lMKISoY3b33A7pw9yXwyRXE8//TR79+5lyZIl3m7F55ggcYDNmzezYsUKwlr0QPwDr7ud0we5L4cIQHx8PHXr1uWdd97xdis+xwSJA7z77rv4B5YjvMV9193G6YPc2+vaFEVAQAB//vOfWbNmDWlpZunrvEyQ+LgjR44we/YcQpr9Dv/QyAK3cfog94V1bYrq8ccfJyIignfffdeSeqWFCRIfN3HiRLKyLlA+rqDpbkvHIPeVdW2Konz58gwZMoTk5GT2799f+BfKCBMkPuzMmTO894/3CWnQpsDJi0rLIPeVdW2K6plnngEwv5XkYYLEh02ZMoUTGb8ReXufaz4rTYPcafXr1q3LgAEDmDx5MmYOYRcTJD7qwoULjBs/npA6zSlXM+aqz5w8CEtDfYAXXniBc+fO8f7779tS32lMkPioWbNmcejgQSLaPXzV+04fhE6vnysmJoYHH3yQf/7zn5yy4IqQ05kg8UEXL17krTFjKFf9ZoKjW1553+mD0On18xs1ahQnT55kwoQJtu/L15kg8UGzZs1iz+7dlL+975U5R5w+CJ1evyBxcXF0796d8ePHc/LkSY/s01eZIPExFy5c4LXXXiO2RStCbr4NcP4gdHr9G3nttdfIyMjgvffe8+h+fY0JEh8zdepU9u3bxzMjRiMijh+ETq9fmNatW9OrVy/eeecdfvvtN4/v31fYGiQicp+IbBeRnSIysoDPB4jIxpw/34jIrXb24+vOnTvHm2++SadOnWh/RxfHD0Kn1y+qV199ldOnTzN+/Hiv9eBttgWJiPgDE4BuQFOgn4g0zbfZHuBOVb0FeB2YbFc/TvDuu+9y8OBB3njjDc6ezXT0IHR6/eKIjY2lX79+/OMf/yizd7vaeUTSFtipqrtVNQuYB1x1n7eqfqOqGTkvv8O1PnCZdOjQId566y3i4+O54447OJi+37GD0On1S+Jvf/sbly9f5sUXX/R2K15hZ5DUBPLGc3rOe9fzOPCZjf34tFdeeYXz588zduxYAKrXqu3IQej0+iVVt25dnn32WRISEsrkk8F2BklBayUUuIiqiPwOV5C8cJ3Ph4pImoiklcZbkjdv3szUqVN58sknadTItVZNWFjRJncuDqcPcl8NkVyjRo0iKiqK5557rsytF2xnkKQDeZ80qwUcyL+RiNwCTAUeUNXjBRVS1cmqGqeqcVFRUbY06y2qyvDhw4mMjOTll1+2bT9OH+S+HiLgWpnv9ddf56uvvipzC2rZGSSpQEMRqSciQUBfYHHeDUSkDvAJMFBVd9jYi89KSEhg9erVjBkzhsqVK9uyD6cPcieESK4nnniCuLg4/vKXv5Spm9RsCxJVvQQMB5YD24BkVd0iIsNEZFjOZi8DlYEPRWS9iJSpk8uMjAyef/552rVrxxNPPGHLPpw+yJ0UIgD+/v5MmjSJI0eOMHr0aG+34zEBdhZX1RQgJd97k/L8+xOAPSPIAUaOHMnx48dZsWIFfn7WZ7rTB7nTQiRX69atefLJJ/nwww8ZNGgQbdq08XZLtjN3tnrJihUrmDx5Mn/5y19o0aKF5fWdPsidGiK53njjDWrUqMEf/vAHzp8/7+12bGeCxAtOnDjB448/TkxMDK+//rrl9Z0+yJ0eIuD64XXq1Kls27aNl156ydvt2M4EiRc8++yzHDx4kJkzZxIcHGxpbacP8tIQIrnuvfdehg4dyt///ne+/vprb7djKxMkHjZnzhxmzpzJqFGjLD93dvogL00hkuvtt98mOjqa/v37l+qH+kyQeND27dv505/+RMeOHXnllVcsre30Qe6EdW1KIiIigqSkJA4ePMjgwYNL7Y1qJkg85OzZszz88MOEhISQmJhIQIB1F8ycPsidtK5NSbRp04a3336bJUuWlNqZ502QeMDly5cZPHgwmzdvJiEhgVq1rHs2sTQMcieta1NSTz/9NL169WLEiBEsX77cq73YwQSJB7z66qvMnz+fcePGcd991192s7hKyyB32ro2JSEizJw5k+bNm9OnTx+2bt3q1X6sZoLEZgkJCbz22ms89thjPPfcc5bVLU2D3In1SyI8PJzFixcTEhLC/fffz6FDh7zdkmVMkNho0aJFPPbYY3Tp0oVJkyZdmcjZXU4fhE6v7446deqwaNEijhw5wt13311qruSYILHJF198QZ8+fYiLi2PhwoUEBQVZUtfpg9Dp9a1w2223sXDhQnbs2EH37t05ffq0t1tymwkSG3z66afcf//9NG7cmJSUFCIiIiyp6/RB6PT6VrrrrrtITk4mLS2Nrl27cvx4gTNoOIYJEoslJiby4IMP0rx5c7788ksqVapkSV2nD0Kn17fDAw88wIIFC9i4cSOdO3fm4MGD3m6pxEyQWERVeeONN+jfvz/t27dn5cqVVKlSxZLaTh+ETq9vp549e5KSksKePXto27YtP/30k7dbKhETJBbIzMykf//+vPTSSwwcOJDly5dT3qLLjU4fhE6v7wldunRh7dq1iAgdOnQgOTnZ2y0VmwkSN/3444+0atWKpKQkxo4da+mDeE4fhE6v70ktWrQgNTWVli1b8sgjjzB06FAyMzO93VaRmSApoaysLN58803atWtHZmYmX375JSNGjLDsEm9m5hlHD0Kn1/eGatWqsWrVKl544QWmTp1K69atWbt2rbfbKhITJCWwevVqWrZsyejRo4mPj2fDhg107tzZ0n2YdW28V9+bgoKCGDNmDJ9//jnnzp2jU6dODB06FF9fPcEESTFs3ryZnj170rlzZzIzM/n0009JTk62ZdJms66Nd+r7iq5du7J161aef/55pk2bRoMGDXjzzTd99nTHBEkhVJVvv/2W+Ph4mjdvzpo1a3jrrbfYunUrPXr0sG2/Zl0bz9f3NWFhYYwfP55NmzbRpUsXRo8eTd26dfnrX//KsWPHvN3eVUyQXEdGRgYTJ06kVatWtG/fnjVr1vDyyy+za9cuRo4cSWhoqLdbLBanD/KyFiJ5xcTEsHDhQr799ls6dOjAq6++Ss2aNenfvz+rVq0iOzvb2y3aO4u80/z6668sW7aMTz75hM8//5yLFy9y6623MmHCBAYNGkR4uPVHCZ7g9EFelkMkr3bt2rFo0SK2bt3KxIkTSUhIIDExkerVq9O7d2969uxJp06dLJ++syjEzhmbROQ+4B+APzBVVcfk+1xyPu8OnAUGq+qPN6oZFxenVqytqqrs2rWL1NRU1q5dy+rVq9myZQvgerCqT58+PPLII7Ru3dqyKzHFsfXAKbq/v8btOk4f5FbVT/7T7bStZ81dxr7i7NmzLF68mPnz55OSksL58+cJCQmhY8eOdOrUiQ4dOtCqVSsqVKhgyf5EZJ2qxhX0mW1HJCLiD0wA7sa1fGeqiCxW1bwTMXQDGub8uQ2YmPNPS6gqJ06cYP/+/ezZs4ddu3axfft2tm7dyqZNm66shBYeHk779u0ZNGgQ3bp1IzY21ivhYTWnDHJv1Xe60NBQ+vbtS9++fcnMzGT16tV89tlnrF69mldeeeXKtI716tUjNjaWmJgYGjZsSP369YmOjqZmzZqUK1fOkl7sPLVpC+xU1d0AIjIPeADIGyQPALPU9V/8nYhUEJHqqnrdhw5Onz7NkiVLyMzMJDMzk9OnT3Py5ElOnjxJRkYGx48f59ixYxw+fJhDhw5x9uzZq75fqVIlmjVrRr9+/WjVqhVxcXE0b97c0qkPfYHTB7kJkeIJCwuje/fudO/eHXD9xpeamsq6detYv349W7duZdmyZVy8ePGq71WpUoVq1apRrVo1KleuTOXKlYmMjKRChQpEREQQHh5OeHh4ob8J2jl6agL787xO59qjjYK2qQlcN0h27NjB73//+2veDwsLo2LFilf+Z7Rr145q1apRq1YtatWqRb169WjQoAEVK1Z0xNFGcKAfsTVLdpv9yX3b2bbuY9o/9Xci61g/xaAT64cG+VtSxykqVqzIPffcwz333HPlvUuXLpGens6uXbvYt28f6enp/Prrrxw+fJgjR46wceNGjh8/zsmTJ68JnMLY9huJiDwM3JuzLCciMhBoq6pP59lmKfCWqq7Neb0SGKGq6/LVGgoMzXnZGLBrJt8qgG9dV7sxp/ULzuvZaf2CfT3XVdWogj6w84gkHaid53Ut4EAJtkFVJwOTrW4wPxFJu96PSb7Iaf2C83p2Wr/gnZ7tvI8kFWgoIvVEJAjoCyzOt81iYJC4tANO3uj3EcMwfJNtRySqeklEhgPLcV3+naaqW0RkWM7nk4AUXJd+d+K6/PuYXf0YhmEfWy9VqGoKrrDI+96kPP+uwFN29lBMtp8+Wcxp/YLzenZav+CFnm29Ic0wjLLBPGtjGIbbTJDkISLjReRnEdkoIgtEpIK3eyqMiDwsIltE5LKI+OzVBRG5T0S2i8hOERnp7X4KIyLTROSIiGz2di9FJSK1RWSViGzL+TvxZ0/t2wTJ1T4HYlX1FmAHMMrL/RTFZqAX8JW3G7mePI9LdAOaAv1EpKl3uyrUDMC69VU94xLwnKrGAO2Apzz1/9kESR6qukJVL+W8/A7XfS0+TVW3qapdN+hZ5crjEqqaBeQ+LuGzVPUrwFHL4KnqwdyHXlX1NLAN153itjNBcn1/BD7zdhOlxPUehTBsIiLRQEvge0/sr3Q9qVYEIvIFcFMBH72oqotytnkR12HiHE/2dj1F6dnHFfRwk7lcaBMRCQf+DTyrqqc8sc8yFySqeteNPheRPwD3A13VR66NF9azAxTpUQjDfSISiCtE5qjqJ57arzm1ySNnIqYXgN+r6tnCtjeKrCiPSxhuypko7CNgm6q+48l9myC52gdABPC5iKwXkUmFfcHbRORBEUkHbgeWishyb/eUX84P2LmPS2wDklV1i3e7ujERSQS+BRqLSLqIPO7tnoqgAzAQ6JLz93e9iHT3xI7Nna2GYbjNHJEYhuE2EySGYbjNBIlhGG4zQWIYhttMkBiG4TYTJIZhuM0EiWEYbjNBYthORNrkzPESLCJhOXNlxHq7L8M65oY0wyNE5A0gGAgB0lX1LS+3ZFjIBInhETnP2KQC54H2qprt5ZYMC5lTG8NTKgHhuJ5lCvZyL4bFzBGJ4REishjXzGj1gOqqOtzLLRkWKnPzkRieJyKDgEuqOjdn/tZvRKSLqn7p7d4Ma5gjEsMw3GZ+IzEMw20mSAzDcJsJEsMw3GaCxDAMt5kgMQzDbSZIDMNwmwkSwzDcZoLEMAy3/X+dsgi17iOQawAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1, figsize=(4, 3))\n", "ax.plot(x, f, color=\"black\")\n", "int_indices = (a < x) & (x < b) # zum \"Herausfischen\" der Indizes, welche innerhalb des Intervalls liegen\n", "ax.fill_between(x[int_indices], f[int_indices], hatch=\"/\")\n", "ax.set_xlim(x[0], x[-1])\n", "ax.set_ylim(0, None)\n", "ax.set_xlabel(\"x\")\n", "ax.set_ylabel(\"f(x)\")\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "ff4b65d4", "metadata": {}, "source": [ "Wir werden am Anfang sehen, wie wir die Integration solcher Integrale explizit durchführen können. Praktisch wird es aber so sein, dass Sie implementierte Routinen verwenden werden, die in Python durch geeignete Bibliotheken zur Verfügung gestellt werden. Nichtsdestotrotz ist es überaus hilfreich zu verstehen, wie diese Routinen funktionieren, um sie sinnvoll anzuwenden.\n", "\n", "Als Zweites sei gesagt, dass die Grundaufgabe immer darin besteht, mit möglichst wenig Rechenaufwand das Integral möglichst präzise zu berechnen. Der Rechenaufwand wird hier verstanden als die Anzahl der Funktionsevaluationen. Für viele der hier besprochenen Funktionen ist der numerische Aufwand nicht besonderns groß und nominell können Sie problemlos sehr viele Funktionsevaluationen durchführen in sehr kurzer Zeit. Dies wird aber in Zukunft nicht immer der Fall sein. Daher ist eine effektive Art der Berechnung besimmter Integrale notwendig." ] }, { "cell_type": "markdown", "id": "0506b56d", "metadata": {}, "source": [ "Ausgangspunkt aller weiterer Betrachtungen ist die Diskretisierung der $x$-Koordinate in $N$ Teilintervalle durch Definition von $N+1$ Stützstellen. Jedes dieser Teilintervalle geht dann von $x_i$ bis $x_{i+1}$ und die Intervallbreite ist gerade definiert als $h = x_{i+1}-x_i$. Dabei ist $x_0 =a$ und $x_N =b$. Jede Stützstelle ist dann definiert als $x_i=x_0+i\\cdot h$. \n", "\n", "![Diskretisierung des Integrationsintervals](Integraldiskretisierung.png \"Diskretisierung des Integrationsintervals\")\n", "\n", "*Diese konstante Diskretisierung vereinfacht unsere Betrachtungen im Folgenden. In der Realität werden Sie tendenziell eher eine adaptive Auflösung verwenden. Bereiche in denen die Funktion sich rasch ändert werden Sie dann feiner diskretisieren. Bereiche in denen die Funktion sich eher wenig ändert werden Sie dann gröber diskretisieren.*\n", "\n", "Das von uns zu berechnende Integral, also die Fläche unter der Kurve, ist dann also zerlegt in die Summe der Flächen in jedem einzelnen Teilintegral\n", "$$\n", "I=\\sum_{i=0}^{N-1}\\int_{x_i}^{x_{i+1}} f(x)\\mathrm{d}x\n", "$$\n", "\n", "Die Aufgabe besteht nun darin, jedes einzelne Teilintegral möglichst präzise zu bestimmen." ] }, { "cell_type": "markdown", "id": "25fb0e5f", "metadata": {}, "source": [ "### Vorbereitung" ] }, { "cell_type": "markdown", "id": "b5e70539", "metadata": {}, "source": [ "Im Folgenden werden wir zwei Funktionen betrachten. Zuerst betrachten wir die einfache Funktion\n", "$$\n", "f(x)=x \\quad .\n", "$$\n", "Das erleichtert uns den Vergleich zum exakten Ergebnis und zur einfachen Visualisierung des Einflusses der Diskretisierung.\n", "\n", "Desweiteren betrachten wir noch das bestimmte Integral\n", "$$\n", "I=\\int_{\\sqrt{\\pi}}^{2\\sqrt{\\pi}}2x\\cos x^2\\mathrm{d}x \\quad ,\n", "$$\n", "Aus der symbolischen Rechnung weiter oben wissen wir, dass dieses bestimmte Integral $0$ ergibt." ] }, { "cell_type": "markdown", "id": "56e437d0", "metadata": {}, "source": [ "Zunächst definieren wir uns also diese beiden Funktionen." ] }, { "cell_type": "code", "execution_count": 30, "id": "195e83f5", "metadata": {}, "outputs": [], "source": [ "def f(x):\n", " return x\n", "\n", "def g(x):\n", " return 2 * x * np.cos(x ** 2)" ] }, { "cell_type": "markdown", "id": "e59da750", "metadata": {}, "source": [ "Dann definieren wir uns die Intervalle, auf denen wir diese beiden Funktionen auswerten möchten sowie die Diskretisierung des Raumes, welche wir für beide Funktionen unterschiedlich wählen werden." ] }, { "cell_type": "code", "execution_count": 31, "id": "812a5c78", "metadata": {}, "outputs": [], "source": [ "# Parameter für f\n", "fa = 0\n", "fb = 2\n", "fN = 8 # Anzahl an Stützstellen für f\n", "fh = (fb - fa) / fN\n", "\n", "# Parameter für g\n", "ga = np.sqrt(np.pi)\n", "gb = 2 * np.sqrt(np.pi)\n", "gN = 50 # Anzahl an Stützstellen für f\n", "gh = (gb - ga) / gN" ] }, { "cell_type": "markdown", "id": "64e51b91", "metadata": {}, "source": [ "Bevor wir nun zur Integration schreiten, werden wir beide Funktionen noch visualisieren. Hierfür wählen wir eine sehr feine Diskretisierung, um eine möglichst glatte Darstellung der Funktion zu erreichen." ] }, { "cell_type": "code", "execution_count": 32, "id": "0eedae50", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAikAAADQCAYAAADCkfHmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAyDklEQVR4nO3deXwV9bnH8c+TjX1NAoQESICA7AECsogr7ijuYG9tq711qVTQtra2Vm977+1ta68K2orc1lqXClgBUVDrLlaQLWHfIiEkISxZSIAQsj33j3OwMQYSyDmZmXOe9+uVF0nOZOaZDDnnOfP7zm9EVTHGGGOMcZsIpwswxhhjjGmINSnGGGOMcSVrUowxxhjjStakGGOMMcaVrEkxxhhjjCtZk2KMMcYYV4pyuoAzFRcXp8nJyU6XYYwB1q1bV6iq8U7XcabsecQY9zjd84jnmpTk5GTWrl3rdBnGGEBEcpyu4WzY84gx7nG65xEb7jHGGGOMK1mTYowxxhhXClqTIiK9RORDEdkmIltEZGYDy4iIzBGRLBHZKCKjglWPMcYYY7wlmJmUauCHqrpeRDoA60TkXVXdWmeZK4FU/8e5wDP+f40xxhgT5oJ2JkVVC1R1vf/zI8A2ILHeYlOBF9RnFdBZRBKCVZMx5vRUldfW5fHYO9udLsUY41H7Syt48r2d/OHDrGavq0Wu7hGRZGAk8Hm9hxKB3Dpf5/m/V1Dv5+8E7gTo3bt30Oo0JpyVHq/i54s38ebGAsamdKWyupaYKIutGWOaRlV5cVUO//XmNqpra5k8qHuz1xn0JkVE2gOvAbNUtaz+ww38iH7tG6rzgHkA6enpX3vcGNM8q7OLuX9BJvvLKvjRZQO458L+REY09OdpjDENe+Ldncz5IIuLBsbzq6lD6dW1bbPXGdQmRUSi8TUoL6vqogYWyQN61fk6CdgXzJqMMf9SVVPLnPd38YcPs0jq0pa/3z2ekb27OF2WMcZjXl2by5wPspiW3ov/uWEYEQF6kxO0JkVEBPgzsE1VHz/FYkuBGSIyH19gtlRVC06xrDEmgHKKjjFzfiaZuYe5cVQSv5w6hPatPDe/ozHGYVkHj/DI61uY0C+WXwewQYHgnkmZCNwGbBKRTP/3fgb0BlDVucBy4CogCygHbg9iPcYYfOPGi9bn88jrm4mIEJ66dSTXjOjpdFnGGA9SVX762iZaR0fw5LS0gA8TB61JUdVPaThzUncZBe4NVg3GmK/6Sjg2uSuPTxtBUpfmjxsbY8LT4ox81uaU8Nsbh9GtY+uAr9/O7RoTJuqGY3946QC+f5GFY40xZ+9EdQ2PvbODEUmduHl0r8Z/4CxYk2JMiLNwrDEmGBauzaOgtILf3TQ8oDmUuqxJMSaEWTjWGBMMJ6preObDLEb17sx5/eOCth17tjImBIVjOFZEOgN/Aobim2/pDlVd6WhRxoSoZRsL2Fdawa9vGIbvYt7gsCbFmBBTeryKh5ds5o0N+8ItHDsbeFtVbxKRGCAsdtoYJ/x1ZQ5949txwYD4oG7HmhRjQki4hmNFpCNwPvAdAFWtBCqdrMmYUJWZe5gNuYf55bVDgnoWBaxJMSYkVPvDsU/7w7Gv3j2eUeEVju0LHAL+IiIjgHXATFU9dnIBuweYMYHxwso9tIuJ5IZR9e8ZHHh29zBjPC6n6Bg3P7uSOR9kcf3IJJbdd164NSjge8M1CnhGVUcCx4Cf1l1AVeeparqqpsfHB/cUtTGhqqyiimUbC7huZCIdWkcHfXt2JsUYjwrHcOxp5AF5qnryTut/p16TYoxpvrc37edEdS03jk5qke1Zk2KMB9UNx45J7sIT09LCJRzbIFXdLyK5IjJQVXcAlwBbna7LmFCzKCOPlLh2jOzVuUW2Z02KMR6zZk8xs+aHXzi2CX4AvOy/smc3di8wYwIq//BxVu0u5oFLBwQ9MHuSNSnGeISFY09PVTOBdKfrMCZULcnIB+D6kcEPzJ5kTYoxHrC3qJyZCzLI2HuYG0Yl8strh7RIaM0YY056Y8M+0vt0oVfXlhtatibFGBdTVRZn5PPI61sQgTm3juTa8A3HGmMcsqfwGNv3H+EXUwa36HatSTHGpSwca4xxi3e27Afg8iHdW3S71qQY40IWjjXGuMnbW/YzLLFTi79RsibFGBexcKwJddmFx9iQe5iconIqa2ro0jaG1O4dGJvclTYxkU6XZxqwv7SCjL2H+dFlA1p829akGOMSFo41oaqiqoa/r8vjhZV72Hng6JffjxCoVd/nraIiuGZET+46vy+p3Ts4VKlpyD+2+oZ6rhjao8W3bU2KMQ6zcKwJZe9s2c9/LN1CQWkFw5M68eg1g5nYP46UuHZERQiHy6vYlF/K21v2s3h9Posz8rljYjIPXDrQzqy4xDtb9tMvvh39u7V882hNijEOKj1exS+WbGaphWNNiCmvrOahRZt4PXMf5/TowO9vHsGEfrFfmwSsS7sYzh8Qz/kD4vnRZQN57J0d/N+KbFbsKuTZ20bTJ7adQ3tgAI6eqGZ1djF3TExxZPvWpBjjEAvHmlCVW1zO915Yy44DR7h/8gC+f1E/oiMbv59t13Yx/M8Nw7h8SHdmzs/k+j9+xovfHcuQnp1aoGrTkH9mFVJVo1w4sJsj27e7IBvTwqprann8HzuY9uxKIiOEV+8ezw8uSbUGxYSErINHuGnuZxSUVvDX28cyc3JqkxqUui4c2I3X751Im+hIbp23io15h4NTrGnURzsO0b5VFOnJzgT4rUkxpgXtLSrn5mdXMueDLK4bmciy+86zq3dMyNi6r4xbnl1FrcLCu8Zz/oD4s15Xclw7Ftw1jo5tornj+TXkFpcHsFLTFKrKRzsOMrF/7Bk3moFiTYoxLUBVWbQ+j6vmrCDr4FHm3DqSx29Js6t3TMjYW1TOt577nFZRESy8azwDezQ/ZJnUpS3P3z6WqhrlO39ZTVlFVQAqNU2188BRCkoruMihoR6wJsWYoCurqGLm/EweWLiBQQkdeGvmJLt6x4SU4mOVfPsvq6mqUV787lhS4gIXdu3frT3P3jaaPUXlPLRoE6oasHWb0/tox0EALhh49mfEmsuaFGOCaM2eYq58cgXLNhXww0sHMP/O8Xb1TpCISKSIZIjIm07XEk4qq2v53gtr2Xf4OH/+dnpQLlMd1zeWH142gGUbC3hldW7A128a9tGOQ5zTowMJndo4VoM1KcYEQd1wbEQEFo5tGTOBbU4XEW7+e9lW1uWU8L+3jCA9uWvQtnP3+f2YlBrHL9/YQnbhsaBtx/gcPVHN2pxiR8+igDUpxgRc/XDs8vsmWTg2yEQkCbga+JPTtYSTJRn5/HVlDt+blMKU4cEdwoyIEH5/8whioiJ4aNFGG/YJss9OXno8wLk8CgSxSRGR50TkoIhsPsXjF4pIqYhk+j8eCVYtxrSUxRkWjnXIk8CDQK3DdYSNXQeO8NCiTYxN7sqDV5zTItvs3rE1P7tqEKt2F7NgjQ37BNNnXxTRJjqS0X2cfYMVzDMpzwNXNLLMClVN83/8Koi1GBNUvnBsBvcvsHBsSxORKcBBVV3XyHJ3ishaEVl76NChFqouNFVW1zJzfiZtYyJ5+hsjW/Ty1OljejGub1d+vXwbJccqW2y74eafWYWMSelKTJSzAy5B27qqfgIUB2v9xrjFWn849s2NBTxw6QBe+d44C8e2rInAtSKyB5gPXCwiL9VfSFXnqWq6qqbHxzs7zu51T763k60FZfzmxuF069i6RbctIvzy2qEcPVHN7Pd3tei2w8XBIxXsOniUCf1inS7F8UzKeBHZICJviciQUy1k74CMG1XX1PL4uzu5pU449r5LUolyaNKjcKWqD6lqkqomA9OBD1T1mw6XFbLW5RQz9+MvuCU9iUsHd3ekhoE9OnDr2N68uCqHrINHG/8Bc0ZWflEEwMR+cQ5X4myTsh7oo6ojgKeAJada0N4BGbfZW1TOLc+uZM77u7guzcKxJjyUV1bzwMIN9Ozchl9MGexoLfdfOoC20ZH8erld0BVo/8wqpGPrKAb37Oh0Kc41KapapqpH/Z8vB6JFxPm2zZhGnAzH7jp4lNnT03h8moVj3UJVP1LVKU7XEapmv7+LnKJyHrtphOP/5+Pat+Kei/rxwfaDrMuxZEEgffZFEeP7xbpiygTHmhQR6SH+e3aLyFh/LUVO1WNMYxoKx05NS3S6LGNaxNZ9ZfxpRTbT0nsx3gVZBYDvTEgmtl0MT7xr2ZRA2VtUTl7JcSb2d8c5g6hgrVhEXgEuBOJEJA94FIgGUNW5wE3APSJSDRwHpqtd+G5cau2eYmbOz2R/WQUPXDqA71/Yz7InJmzU1Co/W7yJzm2ieeiqlrncuCnaxkRx9wX9+O/l21idXczYlOBNJhcu/vlFIYArQrMQxCZFVW9t5PGngaeDtX1jAqG6ppY5H2Tx9Ae7SOzShoV3jXd83gBjWtrLn+eQmXuYJ6el0bltjNPlfMU3x/Xh2U9288S7O3nlznFOl+N5/8wqpFuHVvSLb+90KYDzV/cY41q5xV8Px1qDYsLNwbIKfvf2DialxjE1zX1z/7SJieTuC/qycncR6/eWOF2Op6kqK78oYkK/WPxpDMdZk2JMAxZn5HHl7BXsOmDhWBPefvv2Diqra/nPqUNd88JV361je9OxdRTzPt7tdCmetvPAUYqOVTLBJXkUCOJwjzFeVFZRxS+WbOb1zH2MSe7C47ek0aurTcxmwlNm7mFeW5/H3Rf0IzmundPlnFK7VlF8c1wfnvn4C7ILj5Hi4lrd7PNs37Ur4/u6I48CdibFmC+t3VPMVbO/OnOsNSgmXNXWKv+xdAvxHVox4+L+TpfTqO9MSCY6IoI/rbCzKWdrdXYxPTq2JqlLG6dL+ZI1KSbsVdfU8oR/5lgRWHiXzRxrzOsb8snMPcyDlw+kfSv3n3Tv1rE1149M5O/r8ig8esLpcjxHVVmzx3eFlJuG9exZ2IS1k+HY2RaONeZLx05U85u3tjM8qRM3jkpyupwm+975KZyorrU7JJ+FvcXlHCg7wRiXXcZtTYoJWxaONaZhz3z0BQfKTvDoNUOIcMGso03Vv1sHJvaP5eVVOVTX1DpdjqeszvbN2js22ZoUYxxVd+bYc3p0YLnNHGvMl/JKypm3YjfXpfX05FnF28Yls6+0gg+2H3S6FE9Zs6eYzm2jSe3mjvlRTnL/QKMxAbR2TzGzFmRSUGozxxrTkMff3QnAg1e4Z2bZMzF5UDd6dGzNi6tyuGxID6fL8YzV2cWk9+nqujNn9uxswoKFY41p3Pb9ZSzOyOc7E5Lp2dk9V3iciajICL5xbm9W7Cpk96GjTpfjCQfLKthTVM7YFPedObNnaBPyLBxrTNM89vYO2reK4vsX9nO6lGaZPrYXURHCS6v2Ol2KJ6ze48+jpLhnfpSTrEkxIc3CseFBRHqJyIcisk1EtojITKdr8prV2cW8v/0g91zYz3X35zlT3Tq05vKhPXhtfR4VVTVOl+N6a7KLaRMdyZCeHZ0u5WusSTEhycKxYaca+KGqDgLGAfeKyGCHa/IMVeU3b22je8dW3D4hxelyAmL6mF6UHq/ivW0HnC7F9VbvKWF0ny5Eu3D4230VGdNM9WeOnX+nzRwb6lS1QFXX+z8/AmwDrCtton9sPcD6vYeZNXkAbWIinS4nICb0i6Nnp9YsXJvndCmuVnq8iu37yxjjskuPT7ImxYQMC8caABFJBkYCn9f7/p0islZE1h46dMiR2tyouqaWx97ZQd/4dtw82jsTtzUmMkK4aXQSK3YdYt/h406X41rrcopRhTEuDM2CNSkmROQWlzNt3ioLx4Y5EWkPvAbMUtWyuo+p6jxVTVfV9Pj4eGcKdKHX1ueRdfAoD14+MOQa+ptG90IVFq23symnsjq7hOhIYWQvdz5fhtb/SBOWlmTkc9XsFezcf8TCsWFMRKLxNSgvq+oip+vxgoqqGp54dxdpvTpzeQjOKdI7ti3j+8aycG0etbXqdDmutDq7iGGJnVw7zGdNivGsk+HYWQsyGWjh2LAmvjui/RnYpqqPO12PVzz/2R72l1Xw0yvPcdVN5QLpljFJ7C0u//IyW/MvFVU1bMovdeWlxydZk2I8aV3Ov8Kx90+2cKxhInAbcLGIZPo/rnK6KDcrLa/ijx9mcdHAeMb1de+LVHNdMSSBDq2iWLjWbjpYX8bew1TVqCsncTvJpsU3nlJdU8tTH2Tx1Ae7SOzShoV3jbfsiUFVPwVC81RAkPzx4yyOnKj27PT3TdUmJpKrhyfwxoZ9HL+uxrXDGk5YnV2MCIzu484re8DOpBgPsXCsMYFRUHqc5/+5h+vTEhmU4L4JvAJtaloixypreNfmTPmKNXuKOadHRzq1cW+Gz5oU4wkWjjUmcJ58dxeqcP+lA5wupUWcm9KVhE6teT0j3+lSXKOqppb1e0sYm+zuN3o23GNcrayiikeWbGZJ5j7S+3ThiWlplj0xphmyDh7h1XW5fGdCStj8LUVECNeO6MmfP82m+FglXdt5e9r/QNiyr4zyyhpXh2bBzqQYFzsZjn3DwrHGBMzv3t5B25goZlzc3+lSWtTUtESqa5VlmwqcLsUV1mT7rnZy6yRuJ1mTYlzn5MyxN8/918yxMyfbzLGhQkS6iMgQEekrInZQW9C6nBL+sfUAd53fN+zOJgxK6MCA7u1tyMfv8+xikmPb0q1Da6dLOa0mDfeISDd8l/j1BI4Dm4G1qlobxNpMGMotLmfWgkzW5ZRww8hEfjl1iGVPQoCIdALuBW4FYoBDQGugu4isAv6oqh86WGLIU1V++9Z24tq34o7zQuMmgmdCRJialshj7+wgt7g8rM/K1tYqa3OKuWxwd6dLadRp38WIyEUi8g6wDLgSSAAGAw8Dm0TklyIS+tFw0yIsHBvS/g7kApNUdaCqnuefor4X8Btgqoh819kSQ9sH2w+yek8xMyen0q5VeMYRp6b1BGDphn0OV+KsrENHOVxe5dqbCtbV2P/Uq4Dvqere+g+ISBQwBbgU31TUxpwVC8eGPlW99DSPrQPWtWA5YaemVvnt29tJjm3L9DG9nC7HMUld2jImuQtLMvL5/oX9QnaW3cZ87s+jnOvy0Cw0ciZFVX/cUIPif6xaVZeoaoMNiog8JyIHRWTzKR4XEZkjIlkislFERp15+cbrLBwbXuqfLRGRSBF51Kl6wsWi9XnsPHCUH19+DtFhnu26Ni2RXQePsuPAEadLccya7GK6d2xFr65tnC6lUU363yoiL/rHlE9+nSwi7zfyY88DV5zm8SuBVP/HncAzTanFhAYLx4atS0RkuYgkiMhQYBXQwemiQpnvJoI7GZHUiauGhd5NBM/UlUN7ECGwbGN4XuWjqqzOLmZMcldPnElq6sDkp8DnIvIAkAj8GPjh6X5AVT8RkeTTLDIVeEFVFVglIp1FJEFVw/N/ThipG469fmQiv7JwbNhQ1W+IyDRgE1AO3Kqq/3S4rJD24soc9pVW8PtbRnjiRSnY4tq3Yny/WJZtLOCBSweE3e8kr+Q4+8sqGJvi/jwKNLFJUdVnRWQL8CFQCIxU1f3N3HYiviDdSXn+732tSRGRO/GdbaF3797N3Kxx0pKMfH6xxDcCOHt6mt21OMyISCowE1+ObRBwm4hkqGq5s5WFptLjVTz9YRYXDIhnQr84p8txjauH9eRnizexreAIg3uG17Ufq/15FK80KU0d7rkNeA74Fr5hnOUiMqKZ226ofdWGFlTVef4rAdLj4+ObuVnjhLKKKmbNz2DWgkwG9ujA8pmTrEEJT28Av1DVu4ALgF3AGmdLCl1zP/6CsooqfhLiNxE8U5cP6U5khLBsU/hd5bM6u5hObaIZ0M0bo6xNDQDcCJynqq+o6kPA3fialebIA+rGzJOA8PsfEwZOhmOXbtjHrMmpFo4Nb2NV9X0A9flf4LrmrlRErhCRHf4g/k+bu75QsL+0guc+zea6tMSwO1vQmNj2rRjfN5blm/bjSxyEjzV7ihmT3IWICG8MczWpSVHV61T1YJ2vVwPnNnPbS4Fv+a/yGQeUWh4ltFTX1PLkezu55dlViMCrd49n1uQBFo4NQyJyHoCqltV/TFV3iUhHf5D2bNYdCfwBXxh/MHCriAxuTr2h4Mn3dqIKD4TJTQTP1NXDE8guPMbWgq/9lwxZh46cYHfhMU/Mj3JSY5O5PSwiDe6NqlaKyMUiMuUUP/sKsBIYKCJ5IvJdEblbRO72L7Ic2A1kAf8HfP+s98K4Tm5xOdPmreLJ93Zx7YieLL9vEqP7eOcPwwTcjSLymYg8IiJXi8hYETlfRO4QkReBN4GzvR5yLJClqrtVtRKYjy+YH7Z2HjjCwrW5/Nu43nbW8hQuH9LDN+QTRlf5rNnjrTwKNB6c3QS8ISIVwHr+NZV1KpAGvAf8uqEfVNVbT7di/1U9955hvcYDXs/M5+HFFo41/6Kq94tIF+Am4GagB75bbGwD5jbzCp+GQvjNPdPraf+1bBvtW0Vx38WpTpfiWl3bxTChXyzLNhXw48sHhsVVPquzi2kTHcnQxE6NL+wSjZ13v0lVJwLvAFuASKAMeAnf2PL9qnooyDUajyirqOL+BZnMnJ/JAAvHmnpUtQToCGwE3sU3tUEhcI6IpDVj1U0K4YvInSKyVkTWHjoUuk9bH+44yCc7D3HfJal0CbObCJ6pKcMTyCkqZ8u+8BjyWZ1dzMjenT01oV9jlY4WkT7Av+HLkDwLvIAvje/+qepMizkZjn09M59Zk1NZYOFY07DR+IL3CfhuWHoncCHwfyLy4Fmus0kh/HC4SrCqppb/XraNlLh2fGt8stPluN5lg3sQFSG8GQZDPmUVVWzbX+apoR5ofLhnLvA20BdYW+f7gu+dSt8g1WU8orqmlqc/zOKpD7JI6NSaV+8eb9kTczqxwChVPQrgnxL/78D5+O7f87uzWOcaIFVEUoB8YDrwjcCU6y2vrN5L1sGjzLttNDFR3nm37JQu7WKY2D+O5ZsK+MkVoT3ks25PCaow1kOhWWj83j1zVHUQ8Jyq9q3zkaKq1qCEudzicqbXDcfOtHCsaVRvoLLO11VAH1U9Dpw4mxWqajUwA9+w9DZgoapuaW6hXlNaXsUT7+5kfN9YLh3c3elyPOPq4QnsLS5nU36p06UE1eo9xURFCCN7d3G6lDPS1Bln7wl2IcZbLBxrztLf8N0G43X/19cAr4hIO2Dr2a5UVZfju2IwbD31wS4OH6/i4SmDQvqMQKBdPrgHP4vYxLJNBQxP6ux0OUGzJruYYUmdaBMT6XQpZ8TOB5ozcsTCsaYZVPU/ge8Bh4FS4G5V/ZWqHlPVf3O0OA/LLjzGX1fu4ebRSQzp6Z0rN9ygU9voL4d8QnVit4qqGjbkHfbcUA80/QaDxrAup5hZCzLJLznOrMmpzLiov03MZs6Yqq7Dlz8xAaCqPLp0C62iIvnRZQOdLseTrh6WwIOvbWRTfmlInk3JzD1MVY16LjQLdibFNEHdmWNVbeZYY9zk7c37+WTnIR64dADdOrZ2uhxPumxId6IihGWbQvMqn9XZxYhAugczg/YqY07LwrHGuFd5ZTW/enMrgxI68q3xfZwux7M6t40J6SGfNXuKGdi9A53aRjtdyhmzJsWc0uuZ+Vw1ewU79h9h9vQ0npiWRsfW3vtPbkyomvN+FgWlFfzXdUPszGYzXT0sgdzi42zOD62J3apralmXU+LJoR6wJsU0wMKxxrhf1sEj/GnFbm4enWRnNwMgVId8tuwro7yyxlM3FazLmhTzFetySrhqjs0ca4ybqSoPL9lM25hIfnLlOU6XExI6t41hQggO+azO9t5NBeuyJsUAvlOCs9/bxS3PrrRwrDEu98rqXFbtLuahqwYR176V0+WEjCnDfBO7hdKQz8rdRfSNa0d3j4aq7RXIfBmOfeK9nVwzPMHCsca4WEHpcX69fBvj+8YyfUyvxn/ANFmoDflU19SyOruYcf1inS7lrFmTEubqhmOfnJbGk9NHWjjWGJdSVX6+eDPVtbX85sZhNrNsgIXakM+m/FKOnqhmfF9rUozHNBSOvW6khWONcbOlG/bxwfaD/OiygfSJbed0OSHp6mE92FtczpZ93h/yWbm7CIBx1qQYL7FwrDHec6CsgkeXbiGtV2dun5jidDkh67LBPYgMkSGflV8UkdqtPfEdvJtbsiYljFg41hhvqq1VfvTqBiqqavj9zSOIjLBhnmDp0i6GCf1iWbbR20M+ldW1rN1TwngP51HAmpSwYeFYE4pE5DER2S4iG0VksYh0drqmYPjLZ3tYsauQh68eTP9u7Z0uJ+RNGZ7g+SGfjXmHOV5VwwRrUozbnQzHbrdwrAk97wJDVXU4sBN4yOF6Am77/jJ++/Z2Jg/qxr+d29vpcsJCKAz5rPyiCBE4N8WaFONS9cOxb1k41oQYVf2Hqlb7v1wFJDlZT6CVV1Zz3ysZdGwdxW9uHG5X87SQk0M+Xr7K57MvijinR0e6tItxupRmsSYlRFk41oShO4C3nC4iUFSVhxZtYtfBozwxLc0mbWthVw9LIKfIm0M+FVU1rNtb4ulLj0+yJiXEWDjWhBoReU9ENjfwMbXOMj8HqoGXT7OeO0VkrYisPXToUEuU3iwvrcrh9cx9PDB5AJNS450uJ+xcNsQ35LPcg0M+GXsPU1ld6/nQLECU0wWYwMktLuf+BZmszSnhurSe/Oq6oZY9MZ6nqpNP97iIfBuYAlyipzk3r6rzgHkA6enprj6Hv35vCb96cysXDYzn3ov6O11OWOp68iqfTQX8+PKBnhpqW7nbl0fx6v166rK31yHCwrEmHInIFcBPgGtVtdzpegIht7icO19YS49OrXliWhoRdrmxY7w65LNi1yGGJ3WmUxvvvwZYk+JxRyqqeMAfjk3t3t7CsSbcPA10AN4VkUwRmet0Qc1ReryKO55fQ2V1LX/5zhg6t/V26NHrvDjkU1pexYbcw1yQGud0KQFhwz0eti6nhFkLMsgvOc6syanMuKi/ZU9MWFHVkBkLqaqp5d6X15NdeIwX7hhL/24dnC4p7HWtc5WPV4Z8PvuikFqFSQNCI8dkr2geVFOrzHnfwrHGhIrqmlpmzc/k06xCfn3DMCb0D413waHgqmEJ7CkqZ2uBN4Z8PtlVSPtWUaT16ux0KQER1Fc1EblCRHaISJaI/LSBxy8UkVL/adpMEXkkmPWEAt/MsSt5/F2bOdaYUFBbqzz42kaWbSrg4asHcUt6L6dLMnVc7h/yWbbR/UM+qsonOw8xoV8s0SHypjVowz0iEgn8AbgUyAPWiMhSVd1ab9EVqjolWHWEktcz83l48WYUeHJammVPjPG4mlrlZ4s2sWh9Pj+8dAD/Pqmv0yWZerw05JNdeIz8w8e5+8J+TpcSMMFstcYCWaq6W1UrgfnA1EZ+xjTAwrHGhJ6KqhrufXk9C9bmct/F/ZlxccjEa0KOV4Z8VuwqBOCCEJpXJ5hNSiKQW+frPP/36hsvIhtE5C0RGdLQirw2CVMgrd9bwtVzPmVJZj4zL0ll4V3jbeZYYzyutLyK2/+yhre37OeRKYN54DJ3v0MPd5d75CqfT3Yeok9sW3rHhs5rRDCblIb+4upPoLQe6KOqI4CngCUNrUhV56lquqqmx8eHTod4OifDsTfPXUmtKgvvGs/9l1o41hiv27H/CNf+4VPW5hTzxLQR3HFeitMlmUZ0bRfD+L6xLNvo3nv5VFbXsnJ3EZNC5NLjk4L5ipcH1E2AJQH76i6gqmWqetT/+XIgWkRC6zd8FvJKvh6OTU+2cKwxXqaqLFqfx/V//CfllTW88r1xXD8ypO6HGNKuGeEb8tmQV+p0KQ36PLuI8soaLhzQzelSAiqYTcoaIFVEUkQkBpgOLK27gIj0EP85ThEZ66+nKIg1ud7SDfu4cvYKthUc4YlpI2zmWGNCQNHRE9zz0noeWLiBwQkdefMH59kbD4+5clgCraIiWLw+z+lSGvT+toO0iopgYohdvh60q3tUtVpEZgDvAJHAc6q6RUTu9j8+F7gJuEdEqoHjwPTT3XsjlB2pqOLR17ewKCOfUb07M3v6SMueGONxVTW1vLQqhyff28XxyhoeuvIc/n1SXyJtqnvP6dg6msmDu/PGxgIenjLYVZf4qirvbTvAef3jaBMT6XQ5ARXUGWf9QzjL631vbp3Pn8Y3rXVYW7+3hFnzM8krKWfmJan84GKbOdYYLztRXcOSjHzmfryb7MJjTEqN45Epg0ntbrPIetkNIxNZtrGAj3ccYvLg7k6X86VdB4+SV3Kc718YeleI2bT4DqqpVf7wYRaz399Fj46tWXjXeDsFbIxHqSpb9pXxxoZ9LMrI59CREwzp2ZE/fzudi8/pZlfvhIDzB8QT2y6GxRn5rmpS3tt2AICLzwmtPApYk+KYvJJy7l+QyZo9JUxN68l/XjfUsifGuMi+w8d5YWUOKXFtSY5tR49OrWnXKoqYqAgqqmo4WlFNbslx9hQeI2NvCauzi9lXWkFUhHDBgHhun5jCxP6x1pyEkOjICK4Z0ZO/rd5L6fEq19xl+P1tBxma2JEenVo7XUrAWZPigKUb9vHzxZtQhSemjbCEvzEutKfwGH/+dDdVNY3H5OI7tGJsSldm9IvjyqE96NLO7l4cqq4fmcjzn+1h+aYCbh3b2+lyKDp6gvV7S7jv4lSnSwkKa1Ja0JGKKh5duoVF6y0ca4zbTegfx7ZfXcG+wxVkFx2j8MgJjlVWU1ldS+voSNrGRJLUpS19YtvSrUMrO2MSJoYndaJvfDsWr893RZPy/vaDqMLkQe4Zfgoka1JaiIVjjQkeEfkR8BgQr6qFgVpvVGQEvUNsBk/TPCLCjaOSeOydHeQUHaNPbDtH61m+qYCkLm0YmtjR0TqCxV4lg6zuzLE1tTZzrDGBJiK98N3IdK/TtZjwcOOoJCIjhPlrchtfOIhKy6v4Z1YhVw9LCNkzefZKGUR1Z46dMjyBt2bZzLHGBMETwIN8/bYbxgRFj06tuWhgN15dm0tlda1jdfxj636qapSrhiU4VkOwWZMSJPVnjp1tM8caE3Aici2Qr6obmrBs2N6o1ATeN87tReHRSt73X/7rhJNDPcOTOjlWQ7BZJiXAjp6o5pHXN38Zjn1y2kgbzzamGUTkPaBHAw/9HPgZcFlT1qOq84B5AOnp6XbWxTTLBQO6kdCpNa+syeVKB85klJZX8WlWIbdPTAnZoR6wJiWgLBxrTOCp6uSGvi8iw4AUYIP/SToJWC8iY1V1fwuWaMJQZIQwbUwvZr+/i9zi8ha/UvOdMBjqARvuCQgLxxrT8lR1k6p2U9VkVU3Gd+f1UdagmJZyS3ovBJi/puUz26+tyyMlrh0jQnioB6xJaTYLxxpjTHjq2bkNF5/TjQVrcqmoqmmx7eYWl/N5djE3jkoM6aEesCalWSwca4x7+M+oBGyOFGOa4vaJKRQerWRp5r4W2+Zr6/MQgetHhf5s5daknIWjJ6p5YGEm972SQWq39iy/b5JNbW+MMWFoQr9YzunRgT9/mo1q8PPYtbXKa+vzmNAvlsTObYK+PadZk3KGMvaWcNXsFSzJyOe+S1JZeNd4u3rHGGPClIjw3fNS2HHgCJ9mBf9E3ufZxeQWH+fGMDiLAtakNFlNrfLU+7u4qU449gELxxpjTNi7Nq0nce1b8acV2UHf1kurcujUJporh4b2VT0n2StsE5wMx/6vhWONMcbU0yoqkm+N78PHOw+xZV9p0Lazv7SCt7fsZ9qYXrSJiQzadtzEmpRGWDjWGGNMY749IZkOraOY8/6uoG3jb5/nUKvKN8/tE7RtuI01Kadg4VhjjDFN1alNNHdMTOGdLQeCcjaloqqGv63ey0UDu4VVDtKalAZYONYYY8yZuuO8FDq0jmL2e4E/mzJ/9V4Kj1Zy1/l9A75uN7MmpY764dgFFo41xhjTRJ3aRPPd81L4x9YDZOwtCdh6T1TX8Ownuxmb3JVz+8YGbL1eYK++fnkl5dw6bxX/++5Orh6WwPKZkxhj4VhjjDFn4N8n9aVbh1b8xxtbqa0NzLwpr63Lp6C0gnsv7h+Q9XmJNSn8Kxy7taCMJ6aNYM6tI+nUxsKxxhhjzkz7VlH85Ipz2JB7mMUZ+c1e37ET1Tz53k7SenXm/NS4AFToLWHdpBw9Uc0PF27gvlcy6G/hWGOMMQFw/chERvTqzG/e3k5peVWz1vXsx19w8MgJfjFlUMjfp6chYduknAzHLs7I475LUnnVwrHGGGMCICJC+O/rhlJyrJJHl24+6/XsKTzGvBW7mTI8gdF9wjN+EHZNioVjjTHGBNvQxE7MuLg/SzL3sXxTwRn/fG2t8uO/byA6MoKfXz0oCBV6Q1i9MucfPm7hWGNCjIj8QER2iMgWEfmd0/UYc9K9F/VnRK/OPPj3jew8cOSMfnbuJ1+wZk8Jj14zhIROoX8jwVMJmybljQ37uOLJT74Mx86enmbhWGM8TkQuAqYCw1V1CPB7h0sy5kvRkRHM/eYo2sRE8u9/XcuBsoom/dwH2w/w2Ds7mDI8gRtHJQa5SncLapMiIlf43+FkichPG3hcRGSO//GNIjIq0DWcDMf+oF44NhwDSMaEoHuA36jqCQBVPehwPcZ8RUKnNsy7bTRFR08w7dmV5JWUn3b5T3Ye4p6X1jOkZ0ceu2lE2L9WBa1JEZFI4A/AlcBg4FYRGVxvsSuBVP/HncAzgawhY28JV8+xcKwxIWwAMElEPheRj0VkjNMFGVPfyN5deOG751J0tJJrnvqUtzcXoPrVOVSqa2p55qMvuOP5NfSNb89fbx8bNjcRPJ2oIK57LJClqrsBRGQ+vtOyW+ssMxV4QX1Ha5WIdBaRBFU985RRHTW1yjMfZfHEe7vo0bE1C+4ab9kTYzxKRN4DejTw0M/xPYd1AcYBY4CFItJX678C+NZzJ743Q/Tu3Tt4BRvTgNF9urBkxkRm/C2Du19az9DEjlw2uAfxHVqRW1zOmxsL2FtczlXDevA/Nwy3OIJfMJuURCC3ztd5wLlNWCYR+EqTcqZPLpm5h/n9P3Zy7Yie/Od1Q+1gG+Nhqjr5VI+JyD3AIn9TslpEaoE44FAD65kHzANIT08PzFSgxpyBfvHtWTpjIq+ty+PFVTk8/u5OACIjhPQ+XfjFlMFMHtQt7Id46gpmk9LQb7n+E0NTljnjJ5fRfbrwxozzGJrY0Q62MaFtCXAx8JGIDABigEJHKzLmNKIjI5g+tjfTx/bm2IlqDh+vIrZdDK2jbWinIcFsUvKAXnW+TgL2ncUyZ2VYUqdArMYY427PAc+JyGagEvh2Q0M9xrhRu1ZRtGsVzJdh7wvmb2cNkCoiKUA+MB34Rr1llgIz/HmVc4HS5uZRjDHhQ1UrgW86XYcxJjiC1qSoarWIzADeASKB51R1i4jc7X98LrAcuArIAsqB24NVjzHGGGO8JajnmVR1Ob5GpO735tb5XIF7g1mDMcYYY7wpbGacNcYYY4y3WJNijDHGGFcSrwXhReQQkNPExeMIncsRQ2lfwPbH7Zq6P31UNT7YxQTaGT6PtJRQ+z90OuGyr+Gyn9C8fT3l84jnmpQzISJrVTXd6ToCIZT2BWx/3C7U9scLwul3Hi77Gi77CcHbVxvuMcYYY4wrWZNijDHGGFcK9SZlntMFBFAo7QvY/rhdqO2PF4TT7zxc9jVc9hOCtK8hnUkxxhhjjHeF+pkUY4wxxniU55sUEblCRHaISJaI/LSBx0VE5vgf3ygio5yos6masD8XikipiGT6Px5xos6mEJHnROSg/+ZvDT3utWPT2P545tgAiEgvEflQRLaJyBYRmdnAMp46Rm4XLr/zJu6np/5eTkVEWovIahHZ4N/XXzawjOePKTR5XwN7XFXVsx/47gn0BdAX3y3aNwCD6y1zFfAWIMA44HOn627m/lwIvOl0rU3cn/OBUcDmUzzumWPTxP3xzLHx15sAjPJ/3gHY6eW/Hy98hMvvvIn76am/l9PsqwDt/Z9HA58D40LtmJ7Bvgb0uHr9TMpYIEtVd6vvbqjzgan1lpkKvKA+q4DOIpLQ0oU2UVP2xzNU9ROg+DSLeOnYNGV/PEVVC1R1vf/zI8A2ILHeYp46Rm4XLr/zJu5nSPAfp6P+L6P9H/XDnp4/ptDkfQ0orzcpiUBuna/z+PofQlOWcYum1jref7rtLREZ0jKlBYWXjk1TefLYiEgyMBLfO6O6QvEYuUK4/M5Ps5/g0b+X+kQkUkQygYPAu6oasse0CfsKATyuXm9SpIHv1e/qmrKMWzSl1vX4phAeATwFLAl2UUHkpWPTFJ48NiLSHngNmKWqZfUfbuBHvHyMXCFcfueN7Kcn/14aoqo1qpoGJAFjRWRovUVC5pg2YV8Dely93qTkAb3qfJ0E7DuLZdyi0VpVtezk6TZVXQ5Ei0hcy5UYUF46No3y4rERkWh8LyIvq+qiBhYJqWPkBuHyO29sP73499IYVT0MfARcUe+hkDimdZ1qXwN9XL3epKwBUkUkRURigOnA0nrLLAW+5U9XjwNKVbWgpQttokb3R0R6iIj4Px+L7xgWtXilgeGlY9Morx0bf61/Brap6uOnWCykjpHTwuV33pT99Nrfy6mISLyIdPZ/3gaYDGyvt5jnjyk0bV8DfVyjzrpaF1DVahGZAbyD78qY51R1i4jc7X98LrAcX7I6CygHbneq3sY0cX9uAu4RkWrgODBd/ZFqtxGRV/AlveNEJA94FF/QynPHBpq0P545Nn4TgduATf4xZoCfAb3Bm8fIA8Lld96U/fTa38upJAB/FZFIfC/IC1X1Ta++DjWiKfsa0ONqM84aY4wxxpW8PtxjjDHGmBBlTYoxxhhjXMmaFGOMMca4kjUpxhhjjHEla1KMMcYY40rWpBhjjDHGlaxJMcYYY4wrWZNiWpSIjBGRjSLSWkTaiciWBu79YIwxp2TPI+HDJnMzLU5E/gtoDbQB8lT1fxwuyRjjMfY8Eh6sSTEtzn9fojVABTBBVWscLskY4zH2PBIebLjHOKEr0B7ogO+dkDHGnCl7HgkDdibFtDgRWQrMB1KABFWd4XBJxhiPseeR8ODpuyAb7xGRbwHVqvo3/500PxORi1X1A6drM8Z4gz2PhA87k2KMMcYYV7JMijHGGGNcyZoUY4wxxriSNSnGGGOMcSVrUowxxhjjStakGGOMMcaVrEkxxhhjjCtZk2KMMcYYV7ImxRhjjDGu9P+Wxz9ESEcyFwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Wir definieren zunächst die Anzahl an Punkten auf dem Intervall, die wir für die Darstellung wählen.\n", "# 1000 Punkte werden hier als ausreichend fein betrachtet.\n", "# Beachten Sie, dass dies nicht die räumliche Diskretisierung für die Integration betrifft!\n", "n_fine = 1000\n", "fx = np.linspace(fa, fb, n_fine)\n", "f_fine = f(fx)\n", "\n", "gx = np.linspace(ga, gb, n_fine)\n", "g_fine = g(gx)\n", "\n", "fig, ax = plt.subplots(1, 2, figsize=(9, 3))\n", "ax[0].plot(fx, f_fine)\n", "ax[0].set_xlabel(\"x\")\n", "ax[0].set_ylabel(\"f(x)\")\n", "ax[1].plot(gx, g_fine)\n", "ax[1].set_xlabel(\"x\")\n", "ax[1].set_ylabel(\"g(x)\")\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "a3e663b6", "metadata": {}, "source": [ "### Rechteckregel\n", "\n", "Die einfachste Art besteht nun darin, die Funktion exakt an nur einem Raumpunkt pro Teilinterval zu evaluieren und die Fläche unter dem Teilinterval als ein Rechteck anzunähern.\n", "\n", "![Rechteckregel](Rechteckregel.png \"Diskretisierung des Integrals mit Rechtecken\")\n", "\n", "Beachten Sie hier bitte, dass die Stützstelle in der Mitte des Intervals liegen soll, um das Integral möglichst effizient zu evaluieren. Das können wir uns im Folgenden konkret anschauen. \n", "\n", "Das gesamte Integral ergibt sich dann näherungsweise als\n", "$$\n", "I\\approx=\\sum_{i=0}^{N-1}\\int_{x_i}^{x_{i+1}} f_{i+\\frac{1}{2}}\\mathrm{d}x=\\sum_{i=0}^{N-1}hf_{i+\\frac{1}{2}}=\\frac{b-a}{N}\\sum_{i=0}^{N-1}f_{i+\\frac{1}{2}}\n", "$$" ] }, { "cell_type": "markdown", "id": "1d0c779c", "metadata": {}, "source": [ "Im Folgenden evaluieren wir nun die Rechteckregel für unsere beiden Funktionen.\n", "Wir führen zunächst eine Variable `offset` ein, mit der wir die Stützstelle in jedem Teilintervall verschieben. In der Mitte des Intervalls evaluieren wir die Funktion um $\\frac{h}{2}$ verschoben.\n", "Dann erstellen wir uns eine Liste der Stützstellen und werden die Rechteckregel an diesen aus." ] }, { "cell_type": "markdown", "id": "a775fd5d", "metadata": {}, "source": [ "#### $f(x)$" ] }, { "cell_type": "code", "execution_count": 33, "id": "5ed582fb", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Das Ergebnis unserer approximativen Berechnung beträgt 2.0, das exakte Ergebnis ist 2.\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQwAAADCCAYAAAC8GfzjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAZv0lEQVR4nO3deVTU9f7H8edbFndDxS2XKK/aLXELM69ZLpVLmmWbGWh4jDRLS3Nps+ymmXtpbqk3NTW76k9xKcWuylXRBE3QSFNv5hoKChKgLJ/fH0wdDlEMMMN3Zng/zuHI8P3MzOvb115+5rshxhiUUsoe5awOoJRyH1oYSim7aWEopeymhaGUspsWhlLKbloYSim7eVv1xv7+/iYgIMCqt1eqzIqJiblsjKlVnOdaVhgBAQFER0db9fZKlVkicrq4z9WPJEopu2lhKKXsVmhhiEhDEdkhIvEiclRERhQwRkTkYxE5ISKxItLGOXGVUlayZx9GFjDKGHNQRKoCMSISYYz5Ps+YHkAT21c7YJ7tT6WUByl0hmGMuWCMOWj7/hoQD9TPN6wPsMzk2gf4iUg9h6dVSpVIZmZmiZ5fpH0YIhIAtAb251tUHziT5/FZ/lgqiEiYiESLSPSlS5eKGFUpVRLHjh2jQ4cOJXoNuwtDRKoAa4FXjDEp+RcX8JQ/XDdvjFlojAkyxgTVqlWsw8BKqSIyxjB37lxat27NiRMnSvRadhWGiPiQWxYrjDHrChhyFmiY53ED4HyJkimlSuzChQv07NmTYcOG0bFjR44cOVKi17PnKIkAi4F4Y8yMPxkWDgywHS25B0g2xlwoUTKlVImsXbuW5s2bs2vXLubMmcPXX3/NzTffXKLXtOcoSQcgBIgTke9sP3sDaARgjJkPbAF6AieANCC0RKmUUsWWnJzM8OHDWbZsGUFBQSxfvpzbb7/dIa9daGEYY3ZT8D6KvGMMMMwhiZRSxbZr1y4GDBjAuXPnGD9+PG+99RY+Pj4Oe30901MpD3D9+nXGjBlD586d8fX1Zffu3UyYMMGhZQEWXnymlHKM2NhYgoODiYuLY8iQIUybNo3KlSs75b10hqGUm8rOzmbatGm0bduWhIQENm3axLx585xWFqAzDKXc0unTpxk4cCC7du3iscceY8GCBZTGuU1aGEq5gIBxm+0aZ4zh16P/ISliAWCo2fMVYpp0pe30b50b0EYLQyk3kZ2WTNLWT0g7vpfyDe6g5sMj8fGrW6oZtDCUcgPpp2JI3DKL7PRr+N3/HNXufgwp51XqObQwlHJhOTcyuLLzX6Qe2oyPfyNqPzkB3zq3WZZHC0MpF3X9/DEub55BVtI5qrZ9lOr3DUC8fS3NpIWhlIsxOdkkR31J8p5VeFWpSe1+E6l4S0urYwFaGEq5lMykc1zeNJ0bF45T+c7O1HjgBcpVqGJ1rN9pYSjlAowxpH73FVd2LEa8fPB/ZCyV/97R6lh/oIWhlMUuXLhAwpp3yTgVQ4WA1tTsOQLvqv5WxyqQFoZSFlq3bh1hYWFcT75G9QdeoGqbXuTegsY1aWEoZYHk5GRGjBjB0qVLueuuuyjfejA+NRsW/kSL6cVnSpWyyMhIWrZsyfLly3n77beJiopyi7IALQylSs1v96zo1KkT3t7e7N69m/fee8/h96xwJv1IolQpiIuLIzg4mNjYWMLCwpg+fTpVqrjO4VJ76QxDKSfKyclh+vTpBAUFcfHiRTZu3MiCBQvcsixAZxhKOc3PP//MwIED2blzJ48++igLFy4slXtWOJPOMJRyMGMMn3/+OYGBgURHR7NkyRLWrVvn9mUBWhhKOVRiYiJPP/00ISEhtGjRgtjYWEJDQ1363Iqi0MJQykG2bt1KYGAg69evZ/LkyezcuZNbb73V6lgOpYWhVAmlpaXx0ksv0b17d6pXr87+/fsZO3YsXl6lf4MbZ9OdnkoVYGbEcbvG/XwslhUfjubS2Z+4v+9z9Bw0kp2Xy7PTzucDvPpg0+LGLHVaGEoVQ3Z2FttXzSfi87lUq1GLIR9+RtPW7a2O5XRaGEoV0aVzP7HiwzH8/MNh2nTpzeMvjadilWpWxyoVWhhK2ckYQ9Tm1YQvmIyXjw8hb8ykdaeeVscqVVoYStkhJekSq2e8Sfy3u2japgP9XvsAP/86VscqdVoYShUidncE/571Ftcz0nnsxbfo8MizlCtXNg8wFloYIrIE6AUkGGOaF7C8E7AB+J/tR+uMMe85MKNSlsj4NZX18yfy7dZ1NPjbnTw7bip1GjW2Opal7JlhfAbMAZb9xZj/GmN6OSSRUi7gVFw0K6eM4cqlCzzYfygPPvsi3j7W3uLfFRRaGMaYSBEJKIUsSlnu+vXrvPPOO3wyZQo16jXk5RkrCbijtdWxXIaj9mG0F5HDwHngNWPM0YIGiUgYEAbQqFEjB721Uo5x5MgRgoODOXz4MPf0eIo+Q8ZRvmJlq2O5FEcUxkHgFmNMqoj0BNYDTQoaaIxZCCwECAoKMg54b6VKLCcnh1mzZvH666/j5+dHeHg4Jyo0szqWSyrxrl5jTIoxJtX2/RbAR0Rc8x7pSuXz888/88ADDzBq1Ci6d+9OXFwcvXv3tjqWyypxYYhIXbFduysid9teM7Gkr6uUMxljWLFiBS1atODAgQMsXryY9evXU7t2baujuTR7DquuAjoB/iJyFngH8AEwxswHngCGikgWkA70M8boxw3lspKSkhg6dChffvklHTp0YNmyZdx2m3W/Ed2d2HOU5JlCls8h97CrUi5v27ZthIaGkpCQwKRJkxgzZoxHXobuLGXzdDVV5qSlpfHyyy/TrVs3/Pz8+Pbbb3n99de1LIpIC0N5vOjoaO666y7mzJnDK6+8QnR0NK1b67kVxaGFoTxWVlYW77//Pu3btyc1NZXt27czc+ZMKlasaHU0t6UXnymPdOLECUJCQti3bx/9+/dnzpw5VK9e3epYbk8LQ7mFgHGb7RpnjCH18Fau/OdTpJw3/r1Hs6fh/bT+cG+R3m9E1wLPPSzztDCUx8hOvULi1x+TfvIAFW5pSc2er+JdTc8hdCQtDOUR0o5Hkfj1bExmBtW7hlH1rl6I6C46R9PCUG4t53oaSd8s5Ne47fjWaUzNXqPw9dcLG51FC0O5rYyzR0ncNIOslEtUa/80fh36IV4+VsfyaFoYyu2Y7Eyu7l5Byr61ePvVoU7/yVRocIfVscoELQzlVm5cOs3lTdPJTDhFlRYPUb3LYMqVr2R1rDJDC0O5BWNyuBYdzpVdSylXvhK1+r5NpSbtrI5V5mhhKJd35swZEla/RcbpWCr+rR01u7+MV2U/q2OVSXrcSbksYwwrV64kMDCQ6+ePU6P7y9Tq+5aWhYV0hqFcUlJSEi+++CKrV6+mffv2nA4chE/1elbHKvN0hqFcTkREBIGBgaxdu5aJEycSGRmpZeEitDCUy0hPT2f48OE89NBD3HTTTezfv5833ngDb2+dCLsKLQzlEmJiYmjTpg2zZ89mxIgRvz9WrkULQ1kqKyuLiRMncs8993Dt2jUiIiKYNWuW3rPCRelcT1nm5MmThISEEBUVRb9+/Zg7d67es8LF6QxDlTpjDJ9++iktW7YkPj6elStXsmrVKi0LN6AzDFWqfvnlFwYPHsymTZvo2rUr//rXv2jYsKHVsZSddIahSs2GDRsIDAz8fT/Ftm3btCzcjM4wVLHMjDhu99iMtFQ2zP+A/V+voX7jvzNi0mfk3PI3PvrmhBMTKmfQwlBO9b+jMaz4cAxXEs7Ttd8LdAt5CW8fX6tjqWLSwlBOkZV5g63L5/CfLz+leu2bGTZtObc1D7I6liohLQzlcBdPn2DFh6M5d+J72nV/gj5DXqdCpSpWx1IOoIWhHCYnJ4fdG5azadE0KlSqwqB3P6H5Px6wOpZyIC0M5RBXL11k1bRx/Hgoijvadebpke9Ttbre4t/TaGGoEju4YxNrZ08gOyuLp159n3bdn0BErI6lnKDQwhCRJUAvIMEY07yA5QJ8BPQE0oDnjDEHHR1UuZ60a8msnT2BQzs3E3BHa54Z/SG16t9idSzlRPbMMD4D5gDL/mR5D6CJ7asdMM/2p/Jgxw/uZdW0cVy7kkiP516hy9PP4+WlE1ZPV+gWNsZEikjAXwzpAywzxhhgn4j4iUg9Y8wFR4VUriM9PZ1x48Yx/+OPqd3wNga9O5eGTf8w8VQeyhH/JNQHzuR5fNb2sz8UhoiEAWEAjRrpb6dyNwcPHiQ4OJj4+Hju7RNCr8Gv4Vu+gtWxVClyxLUkBe3dMgUNNMYsNMYEGWOCatWq5YC3VqUhKyuLSZMm0a5dO5KTk9m6dSt9h72lZVEGOaIwzgJ5ryBqAJx3wOsqF3Dy5Enuv/9+3nzzTfr27UtcXBwPPfSQ1bGURRxRGOHAAMl1D5Cs+y/cnzGGxYsX06pVK44ePcqKFSv44osvqFGjhtXRlIXsOay6CugE+IvIWeAdwAfAGDMf2ELuIdUT5B5WDXVWWFU6EhISeP755wkPD6dLly589tlnehm6Auw7SvJMIcsNMMxhiZSlwsPDGTx4MCkpKcycOZPhw4dTrpzeNkXl0gPnCoBr164xcuRIFi1aRKtWrdixYwd33nmn1bGUi9F/OhR79+6lVatWLF68mHHjxrF//34tC1UgLYwy7MaNG7z55pt07NiRnJwcIiMj+eCDD/D11RvcqILpR5Iy6vvvvyc4OJhDhw4xaNAgZs6cSbVq1ayOpVycFoaHCBi32a5xxuRwLWYTV3d9hvhUoNZjb/BNrX/QYtJ/i/R+I7o2KU5M5ea0MMqQrJTLJG6ZRcbp76jYuC01uw/Hq4r+LhBlPy2MMuLX+EiStn6CycmiRreXqNKym96zQhWZFoaHy85IJSliHmnf78L35mb49xqFT/WbrY6l3JQWhgdL/+k7ErfMIvvXK9zUMZib7nkSKedldSzlxrQwPFBO5nWuRi7jWvQGvGs0oG7wNMrX052UquS0MDzMjV9OcnnjdDITf6bqXb3xu38g5Xz0MnTlGFoYHsLkZJOyfy1Xd6/Eq1I1aj/1HhVvbWN1LOVhtDA8wKlTp/hl5etcP/c9lZrdS41uw/CqWNXqWMoD6anhbuy3e1a0bNmSG5dPU7PXKPz7jNWyUE6jMww3lZCQQFhYGBs2bKBz584cvz0E72q1rY6lPJzOMNzQxo0bCQwM5KuvvmL69Ols375dy0KVCi0MN5KamkpYWBiPPPII9erVIyYmhpEjR+oNblSp0b9pbiIqKopWrVqxaNEixo4dy/79+2neXH8fiCpdWhguLjMzk7fffpt7772X7Oxsdu3axeTJkylfvrzV0VQZpDs9XVh8fDwhISHExMQQGhrKrFmz9J4VylI6w3BBOTk5zJ49mzZt2nD69GnWrVvHkiVLtCyU5XSG4WLOnTtHaGgoERERPPzwwyxatIi6detaHUspQGcYLmX16tUEBgayZ88e5s+fz8aNG7UslEvRGYaTzIw4bvfY9NQU1s6ewMEdm2h0e0uGjJ1CWv0AZm3/0YkJlSo6LQyL/XgoilVTx5GSdInuA4bT9ZkX8PLSzaJck/7NtEjmjetsWTKDXes+o1aDAIZ/9AWNmrWwOpZSf0kLwwLnTsazYvJoLp7+kQ6PPEvvwaPxrVDR6lhKFUoLoxTlZGezY81ivl76MZWr+fH8xE/5e9v7rI6llN20MEpJ4oUzrJw6lv8diaFlx248MWIClavpLf6Ve7GrMESkO/AR4AUsMsZMzre8E7AB+J/tR+uMMe85Lqb7MsZwYNs6/m/u+4iUo/+YKdzV9RG9xb9yS4UWhoh4AZ8ADwJngQMiEm6M+T7f0P8aY3o5IaPbSr2axJez3ubI3u00bnk3/Ud/SPXaeot/5b7smWHcDZwwxpwCEJEvgD5A/sJQeRzdt4PVM94k/dcUHgkbx319B+pl6Mrt2fM3uD5wJs/js7af5ddeRA6LyFcicqdD0rmh1NRUXnjhBRaPH0LVGv68OmctnZ4I1bJQHsGeGUZBH7ZNvscHgVuMMaki0hNYD/zhF2GISBgQBtCoUaOiJXUD+/btIyQkhJMnT9L5qcH0GDACb19fq2Mp5TD2/LN3FmiY53ED4HzeAcaYFGNMqu37LYCPiPjnfyFjzEJjTJAxJqhWrVoliO1aMjMzGT9+PB06dCAzM5OdO3fSe/BoLQvlcewpjANAExG5VUR8gX5AeN4BIlJXbLv9ReRu2+smOjqsK/rhhx9o3749//znPwkJCSE2Npb77tNzK5RnKvQjiTEmS0ReAraSe1h1iTHmqIgMsS2fDzwBDBWRLCAd6GeMyf+xxaMYY/jkk08YPXo0lStXZs2aNTz++ONWx1LKqew6D8P2MWNLvp/Nz/P9HGCOY6O5rvPnzxMaGsq2bdvo0aMHixcvpl69elbHUsrpdNd9Ef373/+mefPm7N69m3nz5rF582YtC1VmaGHY6erVq4SEhPDUU0/RpEkTDh06xJAhQ/SMTVWmaGHYYceOHbRo0YJVq1YxYcIE9uzZQ9OmTa2OpVSp08L4CxkZGYwaNYouXbpQsWJFoqKiGD9+PN7ees2eKpv0b/6fOHz4MM8++yxHjx7lxRdfZOrUqVSqVMnqWEpZqswURsC4zXaNMznZpBz4P65Gfo5XxarUfuJdNlcNYvN7O4r0fiO6/uFEV6XcXpkpDHtkJf/C5U0zuH72KJWa/oMa3YbhVekmq2Mp5TK0MMg9CevXI9+QtH0BINR8eCSV7+ysR0CUyqfMF0Z2WjKJW+eQfjyK8g2b4//wSLxvqm11LKVcUpkujLSTB0j86iNyMlLx6zSIam37IOW8rI6llMsqk4WRcyODKzsWk/rdV/jUCqDOU//Et/atVsdSyuWVucK4fv4YlzdNI+vKRard3Re/jiGIt4/VsZRyC2WmMEx2Fsl7V5MctRqvqjWp88wkKjQKtDqWUm6lTBTGsWPHuPj5aG5c/JHKzbtS44EwypWvbHUspdyOR58aboxh7ty5tG7dmqyrF/F/9HX8H35Vy0KpYvLYGcb58+cZNGgQW7dupXv37sQ17o93lRpWx1LKrXnkDGPNmjUEBgYSGRnJ3Llz2bJli5aFUg7gUYWRnJzMgAEDePLJJ2ncuDGHDh1i6NChesamUg7iMYWxc+dOWrRowcqVK3nnnXfYs2cPzZo1szqWUh7F7QsjIyOD1157jS5dulC+fHn27NnDu+++i4+PnluhlKO59U7P2NhYgoODiYuLY+jQoUydOpXKlfUIiFLO4pYzjOzsbKZOnUrbtm25dOkSW7ZsYe7cuVoWSjmZ280wfvrpJwYOHEhkZCR9+/ZlwYIF+Pv/4ZesKaWcwG1mGMYYli5dSosWLTh06BBLly5lzZo1WhZKlSLLZhi/pGQwM+K4XWNTk5NY89G7xO7eym3Ng+g/5kMS6zZg1vYfnZxSKZWXy38kif92F1/MeJO0lKv0GjyaTo+HUs5L71mhlBVctjCup6ex8dMp7N20iroBTQmbuIj6jW+3OpZSZZpLFsbp+MOsnDKGy+dP0+mJQfR47hV8fMtbHUupMs+lCiM7K5OIlfPZvnIeN/nXYeiUpfytZTurYymlbFymMBLOnGLFlDGcORZH0AN9eGzY21SsXNXqWEqpPOwqDBHpDnwEeAGLjDGT8y0X2/KeQBrwnDHmoD2vbYxh78aVhH86BR/fCgx4axat7utRpJVQSpWOQgtDRLyAT4AHgbPAAREJN8Z8n2dYD6CJ7asdMM/2519KSUzgi+lv8EP0f2kWdC/9Rk3ippp1irMeSqlSYM8M427ghDHmFICIfAH0AfIWRh9gmTHGAPtExE9E6hljLvzZi6anpjAlrDeZNzJ4/KXx/KN3f70MXSkXZ09h1AfO5Hl8lj/OHgoaUx/408JIvHCGhk2b8+zYqdRueJudcZVSVpLcScFfDBB5EuhmjBlsexwC3G2MeTnPmM3AB8aY3bbH3wBjjDEx+V4rDAizPWwOHHHUirggf+Cy1SGcSNfPfTUzxhTriII9M4yzQMM8jxsA54sxBmPMQmAhgIhEG2OCipTWjej6uTdPXj8RiS7uc+25+OwA0EREbhURX6AfEJ5vTDgwQHLdAyT/1f4LpZR7KnSGYYzJEpGXgK3kHlZdYow5KiJDbMvnA1vIPaR6gtzDqqHOi6yUsopd52EYY7aQWwp5fzY/z/cGGFbE915YxPHuRtfPvXny+hV73Qrd6amUUr9xmxvoKKWs5/TCEJHuInJMRE6IyLgClouIfGxbHisibZydyZHsWL9OIpIsIt/ZvsZbkbM4RGSJiCSISIGHvz1g2xW2fu687RqKyA4RiReRoyIyooAxRd9+xhinfZG7k/QkcBvgCxwG7sg3pifwFSDAPcB+Z2ayYP06AZuszlrM9bsPaAMc+ZPlbrvt7Fw/d9529YA2tu+rAscd8f+es2cYv59Wboy5Afx2Wnlev59WbozZB/iJSD0n53IUe9bPbRljIoGkvxjiztvOnvVzW8aYC8Z2Aagx5hoQT+7Z13kVefs5uzD+7JTxoo5xVfZmby8ih0XkKxG5s3SilQp33nb2cvttJyIBQGtgf75FRd5+zr4fRkFXk+U/LGPPGFdlT/aDwC3GmFQR6QmsJ/eqXk/gztvOHm6/7USkCrAWeMUYk5J/cQFP+cvt5+wZhsNOK3dRhWY3xqQYY1Jt328BfETEU343gjtvu0K5+7YTER9yy2KFMWZdAUOKvP2cXRieflp5oesnInVtNxhCRO4m9795YqkndQ533naFcudtZ8u9GIg3xsz4k2FF3n5O/UhiPPy0cjvX7wlgqIhkAelAP2PbRe3qRGQVuUcK/EXkLPAO4APuv+3ArvVz220HdABCgDgR+c72szeARlD87adneiql7KZneiql7KaFoZSymxaGUspuWhhKKbtpYSil7KaFoZSymxaGUspuWhhKKbv9PxzOfX6VAMtYAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "offset = fh / 2 \n", "\n", "x = np.arange(fa, fb, fh) + offset # Stützstellen\n", "approx = fh * np.sum(f(x)) # Rechteckregel\n", "\n", "print(f\"Das Ergebnis unserer approximativen Berechnung beträgt {approx}, das exakte Ergebnis ist 2.\")\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(4, 3))\n", "ax.plot(fx, f_fine, color=\"black\") # wir plotten die Funktion wie oben\n", "bars = ax.bar(x, f(x), width=fh) # und die Rechtecke\n", "for bar in bars[::2]:\n", " bar.set_alpha(0.5)\n", "ax.set_xlim([fa, fb])\n", "ax.set_ylim([0, None])\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "f281618f", "metadata": {}, "source": [ "#### $g(x)$" ] }, { "cell_type": "markdown", "id": "6443c893", "metadata": {}, "source": [ "Die Berechnung des Integrals für $g(x)$ erfolgt analog wie oben." ] }, { "cell_type": "code", "execution_count": 34, "id": "6872da14", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Das Ergebnis unserer approximativen Berechnung beträgt 0.25037149932651553, wohingegen das exakte Ergebnis 0 betragen sollte.\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQUAAADCCAYAAABADNcpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAiUElEQVR4nO3deVRV1f//8ee+iIwqDjiLSk5JoiIpKuaQhqZppX0+DpCapVmaDVoOZfb99LXSsrQw06zPLy3nIacCFYe0sBTnCTUVZzEZZBbu/v0h9RUnkHvu3fde9mMt1gI57P066+Cbc/bZZx8hpUTTNO1vJtUBNE2zL7ooaJpWgC4KmqYVoIuCpmkF6KKgaVoBuihomlZAKRWdVqpUSdapU0dF15pWou3ateuKlNL3XtsoKQp16tRh586dKrrWtBJNCHG6sG305YOmaQXooqBpWgG6KGiaVoAuCpqmFaBkoFHTNNu6ePEi//3vf4u0rS4KmubEpJR89tlnjB8/nqysrCL9jL580DQnZTabKdusK6+//jqiZlOqv/BVkX5OnylompN68803SdsXTdnW/8anXThCiCL9nD5T0DQntHTpUj755BPKBHW/r4IAuihomtM5f/48Q4YMISQkhPKdnr+vggC6KGia03n11VfJzs5m3rx5CBfX+/55XRQ0zYls2LCBJUuW8M4771CvXr1itaGLgqY5iWnRR3nu5dcoX6UGpZv14tP18cVqRxcFTXMS+7dv4Ez8AcIiRlCqdOlit2NIURBC+AghlgohjgghDgshWhvRrqZpRSOl5OfvZlC5lj8tHu1pUVtGnSlMB36WUjYCmgKHDWpX07QiiI6O5uKpeB7tOwwXF8umH1k8eUkIURZ4BBgEIKXMAXIsbVfTtKKbNm0aZSv40rzD4xa3ZcSZgj+QCHwrhNgthPhaCOF160ZCiKFCiJ1CiJ2JiYkGdKtpGsCBAweIjo4mtFc4pVyLP5bwNyOKQikgCPhSStkcSAfG3rqRlHK2lDJYShns63vPJeI0TbsPM2bMwMPDg9bd/21Ie0YUhbPAWSnljvyvl3KjSGiaZmVpaWksWLCAvn374lW2vCFtWlwUpJQXgTNCiIb5//QocMjSdjVNK9ySJUtIS0tjyJAhhrVp1N2HkcD3Qoh9QDNgskHtapp2D3PnzqVhw4a0adPGsDYNKQpSyj354wWBUsonpZRJRrSradrd1XhhFtu3b+dStTbUHbfOsHb1jEZNc1Bp+zeAMOEd0MnQdnVR0DQHZDabST+0FQ//Frh4GzPA+DddFDTNAcXGxpJ3LRHPBx8xvG1dFDTNAS1atAhcXPGs18rwtvUajU5MSsmBAwfYvXs3CQkJZGRkUL58efz9/QkNDaVKlSqqI2rFkJeXx5IlS/B4IBiTm6fh7eui4IQuXLjA559/zvfff09CQsL/fUOYQJr/+bJNmzYMHz6cvn37UqqU/lVwFNu2bePChQtUCg63Svv6N8GJpKamMmnSJCIjI7mem8uDDz/Cv/sMo07j5lSsVovIracxZ2dw/a8zZJ3ey5WLsURERPDBBx8wffp0OnfurHoXtCJYtGgRHh4eeNRraZX2dVFwElFRUQwePJiLFy8yZMgQKof+m0rV/W7bzuTmiVv1hrhVb8iwTuPZv30Da76eSpcuXWj31LNsWDiH0hYs0KFZV15eHsuWLaNHjx78XtrDKn3ogUYHl5uby4QJE+jatStXrpemSsQnrK/45B0Lwq2EEASGdmH0V6to92QEv6z4jq5du5KUpOee2asdO3Zw+fJlnnrqKav1oYuCA8vMzKR3795MnjwZ76ZhVH12Gm7VGtx3O6Xd3Hnqpbfp/+YUtm3bRmhoKJcuXbJCYs1Sq1atolSpUnTr1s1qfeii4KCSkpLo0qULq1ev5osvvqBi15GYXN0sajO4cy8q9H6Pw8dOUKtJCLVGLTQorWaUH3/8kfbt2+Pj42O1PnRRcECTV+4isGU7Yn//nYgJn5LToIthbbvXDsT36Xe4fvUcl5dOKvJLSTXri4+P58iRI/Tq1cuq/eii4GDS09P5+u2hXDgVz+CJX9DsEeNPIz3qNMP3iTHknD/K0KFDkVIa3od2fz5dH8/oj+cCcNknoNjLtxeFYUVBCOGSvxzbGqPa1ArKzc2ld+/enDq8h/CxH9O4VQer9eXZsA3lQgcwb948PvvsM6v1oxXdwd9iqObfkApVa1q1HyPPFEahV3G2qlGjRhEVFcUzo/6Hpo90tXp/5dr8m169evHWW2+xe/duq/en3V1aylVOHorjodaPWr0vo977UBPoDnxtRHva7SIjI5k5cyajR48mpNszNulTCBNz587F19eX/v37k5GRYZN+tdsd/n0r0mwmIKSj1fsy6kzhM+BNwFzIdloxVO33ASNGvoJHvVYsFu1s2vd3cX/Rc9Rkjhw5Qlj/YVa9ltXu7sjOX/D2qUjN+g9ZvS+Li4IQogdwWUq5q5Dt9BLvxXDp0iWurJ5KqfLVqdTjDYTJxeYZGjRvTbsnI9i+6ntOHdKXEbaWl5dH/K5tNAoOxWSy/r0BI3poC/QUQpwCFgKdhBDzb91IL/F+//Ly8ggPD8ecnYHvk2Ot8kRcUXUb9CrlKlZh8WcTuX79urIcJVFcXBzpqck0bBFqk/6MWM15nJSyppSyDtAXiJFSWufxrRJm8uTJbNiwgfKdh1Hat47SLO6e3jw94h0unopn2rRpSrOUND///DNCCMcpCpp1bNu2jUmTJhEeHo53oHGTkyzxUJvOBLTuxPvvv6+nQdtQVFQUNesH4O1TwSb9GVoUpJSbpZQ9jGyzJEpLS2PgwIHUrl2bmTNnIoRQHekfT7zwJllZWUycOFF1lBIhKSmJ3377jUbBthtg1mcKdubT9fGE9RvKyZMnefzl//B17AXVkQqoXLMuL7/8Ml9//TUHDhxQHcfpbdy4EbPZTMNg21w6gC4Kdufozm38umYBjzw9iAcCH1Yd546WyRBw9eDhnoOoM3at6jhOLSoqinLlylH7wWY261MXBTuSnJzMomkTqFzLn26DXlUd565cPMpQtvW/yDq5i6yz+g2B1iKl5Oeff6Zz5864uNhuPSRdFOzI2LFjSbl6mX5jPqK0m7vqOPdUpnl3TJ4+pGz7XnUUp3Xo0CHOnj1LWFiYTfvVRcFObN++na+++opHnnyW2o0CVccplKm0O+Va9Sbr9F5++eUX1XGcUlRUFIAuCiVRTk4OQ4cOxc/Pj64DX1Edp8i8m3fDxas87777ruooTmnjxo00bNgQP7/Cl9Yzki4KdmDq1KkcOnSIyMhI3Dy8VMcpMpOrO2VD+rBp0ya2bt2qOo5TqT3mR37asIlzng/YfDBXFwXFagydzdvvvodnw7aM2GY/8xGKyrtpV3x9ffnoo49UR3EqORePI3Mycfez/aWkLgoKSSm5Gj0T4eJK+UeHqo5TLCZXN0aOHMm6des4ePCg6jhOIythHwDufk1s3rcuCgotWbKErNN7Kd9+IKXKVFQdp9heeuklPD09+fjjj1VHcRpZp/fi6lsHF89yNu9bFwVFMjIyGD16NK6V6+LdzPqrKFnTd3F/EdTlaebN/553F+g7EZbKysoi+9xh3Gs3VdK/LgqKTJ06lTNnzlCh8zAlayQYrf3TgzCb8/hl5Xeqozi82NhYZG4O7rXV3JrWRUGBhIQEPvroI/71r3/hXsv6K+nYQsVqtWjaLozf1izk2rVrquM4tJiYGBAmZb8buigo8NZbbyGlZMqUKaqjGOqRpwaSlZHG/Pm3rbGj3YeYmBhKV62HyU3N7WkjlmOrJYTYJIQ4LIQ4KIQYZUQwZ/Tp+nhGTvuBhQsX8kifISyPz1YdyVC1H2xGzXoBREZG6ndFFFNaWho7duxQdukAxpwp5AJvSCkfBEKAl4UQjQ1o1+mY8/JY8eX/4lOpKh2feV51HMMJIQjtNYCDBw+yZcsW1XEc0rZt28jNzcXdT80gIxizHNsFKWVc/ufXuPHuhxqWtuuMfo9ezrnjh+jxwhjcPNStt2hNzTp0p0KFCkRGRqqO4pBiYmJwdXXFreaDyjIYOqYghKgDNAd23OF7JXo155SUFNZ9+yl1A4Jo3qG76jhWU9rNnSFDhrBixQrOnj2rOo7D2bRpE61bt8bkqu4pWSNfG+cNLANelVKm3vr9kr6a83/+8x/SU67y5PAJdrW8mjUMHz4cs9nM7NmzVUdxKElJScTFxdGpUyelOYx6Q5QrNwrC91LK5Ua06UyOHj3K9OnTaRnWm1oNnOMW5L3UrVuX7t27M3v2bL0c/H3YunUrZrOZjh2t/xaoezHi7oMA5gKHpZR67e87eOONN/Dw8ODxwa+pjmIzw4YN49KlS6xdq5drK6qYmBg8PDxo1aqV0hxGrPHUFogA9gsh9uT/23gp5ToD2nZ4P/30E2vXrmXq1Km4lK+kOo5N1Bm7Fmk24+Jdgf6jJ1M51pVTHzrvOIpRYmJiCA0Nxc3NTWkOi4uClHIb4NwXycVUe8yPnP9mBKXKV+fzS/V4VXUgGxImF7yadCY1dim5166ojmP3Ll26xIEDBxgwYIDqKHpGozVdi1tL7tWzlO/0PMLFVXUcm/Nu0gWkmfT9G1VHsWufro9nzIwFAFz0qqf8Jb66KFhJYmIiydt/wL1uCzwesM+l2q3NtXw13GsHkrYvGrNZv5D8Xo7ticXd05ua9dXP+9NFwUrefvtt5PUsKnR63ulvQd6Ld+Bj5KZcYtOmTaqj2LXje2J5IPBhmy7lfje6KFjBnj17mDNnDmWCeuBaqZbqOEp5NmiDyd2buXPnqo5it5IuX+DK+dPUaxaiOgqgi4LhpJSMGjWKihUrUq5tP9VxlBOlSuMV0JFly5bx119/qY5jl47vvTEBuL4uCs5p6dKlbN26lffffx8Xd2/VceyCd5Mu5OTksGjRItVR7NLxPbF4lStP1ToNVEcBdFEwVGZmJqNHj6Zp06Y8/7zzPQVZXKWr+BMYGMh33+lVmW4lpeTYnljqBbbEZLKP/472kcJJPPXiWyQkJNA2YjQzYk6ojmNXnn32WXbs2MHRo0dVR7ErJ06cIDnxgt2MJ4AuCoZJSEggZtEcmj7SlXqBLVXHsTv9+/fHZDIxb9481VHsSkxMDGA/4wmgi4JhxowZA8ATL7ypOIl9qlatGo899hjz5s3TcxZuEhMTQ9mKlfGtWVd1lH/oomCALVu2sHjxYjr96wUqVNHry9zJp+vjqdS8CwkJCYycpt9UDTfGE2JiYqjfLMSu5rLoomCh3NxcXnnlFfz8/Oj4zBDVcezaQ60fxc3Ti50bVqqOYhcOHjxIYmKiXV06gDFPSZZoc+bMYd++fSxZsoQz7h6q49i10u4eNG3Xlb1bfyIjIwNPT+dckq6o/h5PsKdBRjBukZWuQoijQojjQoixRrTpCK5evcrbb79Nhw4d6N27t+o4DiG4y5NkZ2awcuVK1VGUi4mJoW7dunZ3yWnEIisuQCTQDWgM9CspqzlPnDiR5ORkpk+fblfXhPbM/6FgylepUeLnLOTl5bF582YeffRR1VFuY8SZQkvguJTyTyllDrAQ6GVAu3at+pBIImd+iVezbvT84Qx1xuoVhorCZDLRotMTrF+/nvPnz6uOo0xcXBwpKSnK12O8EyOKQg3gzE1fn8XJl3iXUnJ1/SxMbl6UCw1XHcfhBHfuhdlsZsGCBaqjKPP3eII9FgVh6Zt8hBDPAGFSyufzv44AWkopR96y3VBgKED5ytVbvDP//x6lnb7x2G3tjnq0foGv7WmbnRtW8sOUt5g1axbDhg27bVutcG7VGyDz8qg+eMY//2bPx7y429z6fYDXujQgLCyMc+fOceDAgdu+b01CiF1SyuB7bWPEmcJZ4Obng2sCt50X3rzEu1e58gZ0q0bGtRRWffURrVq14oUXXlAdx2F5BXTi+uU/yUk8pTqKzeXk5PDLL7/Y5VkCGFMU/gDqCyHqCiFKA32BVQa0a5fWfvMJ6deSmTVrlt08wOKIvBq1A2Ei/eBm1VFsLjY2lszMTLscZARjXhuXC4wAorjxyrjFUsqDlrZrj04f3kvsusW06xVBs2bNVMdxaC5ePnj4tyD90GakLFnTnmNiYjCZTLRv3151lDsy5E+dlHKdlLKBlPIBKeX/GtGmvcnLy2XpjHcpU8GXrgNfUR3HKXg17kDetStkJ9j2ulq1mJgYgoKC8PHxUR3ljvT5bxG5Hl3PuROH+frLLxjXK0h1HId36sPunF44iTJlytDd83iJeS9EdmYGsbGxdjueALooFEluaiLvvPMOYWFh9OnTR3Ucp+Hp6Unv3r1ZunQpmZmZquPYxMmDu7h+/brdjieALgqFklLyV9QX5OXlMXPmTD1z0WDh4eGkpqayevVq1VFs4tieWFxdXWnbtq3qKHeli0Ih0g/GkPXnLiZPnoy/v7/qOE6nQ4cOVK9enfnz56uOYhPH98QSEhKCl5eX6ih3pYvCPeSlJZG0cQ5uNR5kxIgRquM4JRcXFwYMGMBPP/1EWspV1XGsKuNaCmePHbTrSwfQReGerq7/EvP1bCp2ewUXFxfVcZxWeHg4ubm57Nnyk+ooVnVi3x9IKe16kBF0Ubir9CPbyIj/FZ/Q/rhWLNkvdLG2wMBAmjRpwq6NTjvnDcgfT3BzV/6q+cLoonAHv77SHPMvswkODiZxw9cl5naZShEREZw+vIfEc6dUR7GaY7t/xf+hYEqXLq06yj3ponALs9nM4MGDyczMZP78+ZQqpRensoV+/fohhGDXRue8C5F0+QKXEk7QsIX93nX4my4Kt9i+6nuio6OZNm0aDRs2VB2nxKhZsyb1mrZiV8wqLH1y1x7Fx20HoGGLUMVJCqeLwk0unjrG6q+n0r17d/1ItAItOvfir/MJnD68R3UUwx3dtZ2yFXypWuf2R6ntjS4K+a7nZDP/w9G4e3gxd+5cPUlJgcC2j1GqtJvTDTia8/KIj/uVhi1CHeL3SheFfCsi3+f8n0foO+ZDqlSpojpOieTu5c1DrR9lz5Z15F7PUR3HEK91acAjFdPIuJbMq4P68FoX+3iJ7L1YVBSEEFOFEEeEEPuEECuEED4G5bKpCud/I/anxYwfP545E/TCKaq81qUB7495ifTUZBrm/qk6jmGio6MB6NKli+IkRWPpmcJ64CEpZSAQD4yzPJJt5SSeYvjw4XTs2JH33ntPdZwS77HHHsPX19eppj1HR0cTFBSEr6+v6ihFYlFRkFJG5y+yAhDLjaXYHIY5K43ElR/g4+PDDz/8oG8/2gFXV1f69u3LqlWrMGenq45jsWvXrvHrr7/y2GOPqY5SZEaOKTwHOMw8VWnOI/HHj8hNvsSiRYuoWrWq6khavvDwcLKzs0k/sl11FItt2rSJ3NxcwsLCVEcpskKLghBigxDiwB0+et20zQQgF7jrm0OFEEOFEDuFEDvTU5KMSW+BpI2zyTq1m4phL9GuXTvVcbSbPPzww9SvX5/0Q5sK39jORUdH4+XlRevWrVVHKbJCi4KUsrOU8qE7fPwIIIQYCPQABsh7zDqxp9Wcr8Wt4VrcWsq2fBrvQMc5rSsphBBERESQnbCf3NTLquNYJCoqivbt2+Pm5qY6SpFZevehK/AW0FNKmWFMJOua1jKb5I2zeeKJJ7j662L9XIOdGjBgAADph7YoTlJ816+e4/jx4zz++OOqo9wXS8cUvgDKAOuFEHuEELMMyGQ18XG/0q9fP1q1asWCBQv049B2zN/fH7caD5J+YJPDTnvOPPEHAN27O9YfHouG26WU9YwKYm2nD+/lm0kv06BBA9asWWPXK99oN3gFdORq9EyuX/6T0lUeUB3nvmWe+J2AgADq1KmjOsp9KREzGk8ejOOrcc9RpnxFoqOjqVChgupIWhF4NgoFUynSDjregKM5O4OsMwfp0aOH6ij3zemLQlDpC3zz9vPUrF6V3Tu2U61aNdWRtCJy8SiLxwPBZBzagjTnqY5zXzJPxoE5TxcFe5NxfAfdunXDz8+PrVu3UquWXkHJ0XgFdCQvPYms03tVR7kvmSd2YnL3JiQkRHWU++aURUFKSeofP5K47H0aN27M5s2b9RmCAzr1YXfOLJhIuXLlCHM75jB3iqQ0k/nnH7j7t3DIWbJOVxRk3nWurv+SpJg5eNRvxZYtW6hcubLqWFoxubu788wzz7B8+XLS0x1j2nPOhWOYM1LwfKCl6ijF4lRFIWbog1Te/L+k7V7HmDFjSDuyXd9lcAIRERGkp6ezcuVK1VGKJPP47yBMuNd1zNcLOkVRkFLyR/QKmjdvzvHjx1m+fDlTpkzRr4p3EqGhofj5+TnEk5OnPuxOzdQDtH+kHWem91Udp1gc/n/NlfMJfDX2ORZ8PJYmTZqwe/dunnrqKdWxNAOZTCbCw8OJjo4m9Wqi6jj3dPToUQ4ePEjv3r1VRyk2hy0KeRkpJKz7ko+HdufC8f1ERkaydetWh5soohXNgAEDMJvN7N68VnWUe1q2bBmAQ/9hcrih0etJF7gWt4a0fdHMyM1m8ODBvPfee9SoUUN1NM2KGjduTFBQELs2rqL904NUx7mr5cuX06pVK2rWdKilRQpwiKKQlnyVa3t+JuPodrJO7QGTCc+GoexaPotGjRqpjqfZSEREBK+99hoXTsZTra79rXV49eJZdu3axZQpU1RHsYiSonAt6Qr7tkVToWpNPMuUIy8zFSFMmHNzkDmZ5KZc5re1cZw7cYiTh3Zz8WQ8Ukr8/f0Jn/gOw4YNo3r16iqiawqFh4cz+s03+W3dIp5++R3VcW6zb9t6AJ5++mnFSSyjpCikXLnEf/9n5D23WQKULVuWkJAQXhzYn549e9K0aVOHWCJbs45KlSrRtF1Xdq5fSffn3sDNw1N1pAL2b4+madOmPPCA4z28dTMlRaFZs2bMnTuXkydPkpKSQlpaGrm5uXh5eeHl5YWfnx9169alRo0a+raiVkCbHn2Ji1nN7s1rCen2jOo4/0hOvMipQ7uZNGmS6igWM6QoCCFGA1MBXynllcK2d3FxISgoiKAgx5zcoalTN6AFVes04Le1C+2qKMRtWoOUkv79+6uOYjGL/wwLIWoBXYAEy+No2r29/lhD3h79CmfiD9CufKrqOP/YFbMKv0ZNqVfPYZYYuSsjzs0/Bd4EHHN5HM3hhIeH4+npyaxZ9rHQ1/mTR7nw51FaPNpTdRRDWLpGY0/gnJSy0Odab17NOTHRvmelafatXLly9O/fnx9++IG8rDTVcYiLWY3J5ELz9o61FuPdWLrE+wRgYlE6unk1Z0d5U45mv0aOHElmZiZpe35WmsNsNhMXs4ZGD7fD28c5VvQqdKBRStn5Tv8uhGgC1AX25t8mrAnECSFaSikvGppS024RGBhI586d2RS7irIP90K4uCrJEVT6AsmJF5g5/RP6OcDLY4ui2JcPUsr9UsrKUso6Uso6wFkgSBcEzVZef/118tKukn5km7IMs2fPxsfHhyeffFJZBqPpSQCawwoLC8O1Yi2u/bFSyTLweRkpLF++nGeffRYPDw+b928thhWF/DOGQucoaJpRTCYTZYJ7kXPpBNkJ+23ef9r+jeTk5DB06FCb921N+kxBc2heAR0xefqQErvEpv1KKUnb+zNt27YlICDApn1bmy4KmkMzubpRtuVTZJ3aTfa5wzbrN+v0XnKTzjNs2DCb9WkruihoDu3Uh905u+bzGw9LJW6w2YrPqX+swOTlwzPP2M9Ua6PooqA5PC8vL0aPHk1UVBQ7duywen85iafI+nMXZYOewN3d3er92ZouCppTePnll6lYsSITJxZpLp1FUv9YiXB1w7t5N6v3pYIuCppT8Pb2Zty4cURHR3Nk5y9W6yc39TLpBzfj3aQzLh5lrdaPSrooaE5jxIgR+Pv7s2r2FMx5xr978tSH3el2/VdKlzJxYMVMh3lj1f3SRUFzGm5ubnz44YdcPBXPjqhlhrd/7Ngxvv32W1588UX8/PwMb99e6KKgOZU+ffpQNyCIdd98QlryVUPbnjhxIm5ubowfP97Qdu2NLgqaUxFC0OeV98hMT+PHrz4wrN0T+/5g4cKFvP7661SpUsWwdu2RLgqa05kytAcTxo9j18ZVNMo9YXF7ebnXWfb5e9SuXZtx48YZkNC+6aKgOaUJEyYQEBDAoEGDyE2z7DIiZtEcLp4+xowZM/D0tK8VpK1BFwXNKbm5ubF48WLS0tK4snoq0ly8uxGnDu0mat4XNO/QnZ49nWO5tcIYsXDrSCHEUSHEQSGEY78aR3MqjRs3ZubMmWQn7Ofqhq/u+/Hqa0lXmPfBG/j4VqXPqPeslNL+WLTEuxCiI9ALCJRSZgshKhsTS9OMMXDgQF75cjWpO5bh4lGWcqEDivRCIXNOJqunTiL7WhJbtmwhOLiFDdLaB0vPFIYDH0opswGklJctj6Rpxkr+bQnPPfccKb8upGdGFMffD7vn9nkZKVxaOJ64uDgWLlxIcHCwjZLaB0uLQgOgnRBihxBiixDiYSNCaZqRhBDMmTOHN954gy+++IK2bdty7sTtj1lLKck4/jsXvh3J9cTTrFixgieeeEJBYrUKvXwQQmwAqt7hWxPyf748EAI8DCwWQvjLO1y8CSGGAkMBp54Nptknk8nExx9/TMuWLRk+fDi/D3+S+s1bU79Za7zKlSf590NknvidnEsncK3kh2/viSVmYPFWwpK17YQQP3Pj8mFz/tcngBAp5T1f7BAcHCx37txZ7H41zRLJyclERkYyf/58jhw5Atw4m2jRogXPP/88gwYNws3NTXFK6xBC7JJS3vN6yNKi8CJQXUo5UQjRANgI+N3pTOFmuiho9iI1NZXk5GQqVapUIuYgFKUoWPqC2W+Ab4QQB4AcYGBhBUHT7EnZsmUpW9Y5H4EuLouKgpQyBwg3KIumaXZAz2jUNK0AXRQ0TStAFwVN0wrQRUHTtAIsuiVZ7E6FSARO27jbSkBJfa2d3veS5277XVtK6XuvH1RSFFQQQuws7P6ss9L7XvL23ZL91pcPmqYVoIuCpmkFlKSiMFt1AIX0vpc8xd7vEjOmoGla0ZSkMwVN04rAaYqCEKKWEGKTEOJw/nqRo+6wjRBCzBBCHBdC7BNCBKnIarQi7nsHIUSKEGJP/of138RqA0IIdyHE70KIvfn7fttiik583Iuy7/d/3KWUTvEBVAOC8j8vA8QDjW/Z5nHgJ0BwY2GYHapz23DfOwBrVGe1wr4LwDv/c1dgBzfW9CgJx70o+37fx91pzhSklBeklHH5n18DDgM1btmsF/CdvCEW8BFCVLNxVMMVcd+dUv6xTMv/0jX/49aBMmc97kXZ9/vmNEXhZkKIOkBzblTOm9UAztz09Vmc7D/PPfYdoHX+qeZPQogA2yazHiGEixBiD3AZWC+lLDHHvQj7Dvd53J2uKAghvIFlwKtSytRbv32HH3Ga2y+F7HscN6a4NgU+B1baOJ7VSCnzpJTNgJpASyHEQ7ds4rTHvQj7ft/H3amKghDClRv/Kb6XUi6/wyZngVo3fV0TOG+LbNZW2L5LKVP/PtWUUq4DXIUQlWwc06qklMnAZqDrLd9y2uP+t7vte3GOu9MUBXHjDR9zgcNSyml32WwV8Gz+aHQIkCKlvGCzkFZSlH0XQlTN3w4hREtuHPu/bJfSOoQQvkIIn/zPPYDOwJFbNnPW417ovhfnuFu6RqM9aQtEAPvzr7EAxgN+AFLKWcA6boxEHwcygMG2j2kVRdn3PsBwIUQukAn0lfnD0w6uGvD/hBAu3PiFXyylXJO/qLCzH/ei7Pt9H3c9o1HTtAKc5vJB0zRj6KKgaVoBuihomlaALgqaphWgi4KmaQXooqBpWgG6KGiaVoAuCpqmFfD/Af7gna6jQUxRAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "offset = gh / 2 \n", "\n", "x = np.arange(ga, gb, gh) + offset\n", "approx = gh * np.sum(g(x))\n", "\n", "print(f\"Das Ergebnis unserer approximativen Berechnung beträgt {approx}, wohingegen das exakte Ergebnis 0 betragen sollte.\")\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(4, 3))\n", "ax.plot(gx, g_fine, color=\"black\")\n", "bars = ax.bar(x, g(x), width=gh)\n", "for bar in bars[::2]:\n", " bar.set_alpha(0.5)\n", "ax.set_xlim([ga, gb])\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "61d801fe", "metadata": {}, "source": [ "Es gibt auch andere Arten der Diskretisierung. Sie können, z.B., die Fläche in jedem Teilinterval mit einem Trapez annähern. Oder Sie können auch Polynome höheren Grades und nicht nur Gerade verwenden, um die Fläche zu berechnen. Das erfolgt im Rahmen der Simpsonregel. Glücklicherweise müssen Sie sich um diese Details nicht selbst kümmern, sondenr können vordefinierte Routinen verwenden. " ] }, { "cell_type": "code", "execution_count": 35, "id": "62d85faf", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Das Ergebnis unserer numerischen Berechnung mit der vorhandenen Funktion simpson beträgt -0.0005916662694904901.\n" ] } ], "source": [ "from scipy.integrate import simpson\n", "\n", "x = np.linspace(ga, gb, gN)\n", "approx = simpson(g(x), x)\n", "\n", "print(f\"Das Ergebnis unserer numerischen Berechnung mit der vorhandenen Funktion simpson beträgt {approx}.\")" ] }, { "cell_type": "markdown", "id": "5bdfb711", "metadata": {}, "source": [ "Eine sehr häufig verwendete Routine ist die Gauss-Quadratur. Ausgangspunkt der Gausschen-Quadratur ist die Idee, eine Funktion $f(x)$ nicht näherungsweise an einer endlichen Anzahl von Stützstellen zu evaluieren, sondern die Funktion im gesamten zu integrierenden Interval als eine Summe von $n$ Basisfunktionen $P_i(x)$ zu entwickeln, in denen jede einzelne Basisfunktion mit einer bestimmten Amplitude $a_i$ beiträgt! Numerisch brauchen Sie sich über diese Details keine Gedanken machen. Sie können auch hier wieder vorhandene Funktionen nutzen, um das Verfahren effizient zu nutzen." ] }, { "cell_type": "code", "execution_count": 36, "id": "1cb9d668", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Das Ergebnis unserer numerischen Berechnung mit der vorhandenen Funktion quad beträgt 1.7710371678247452e-15.\n" ] } ], "source": [ "from scipy.integrate import quad # siehe https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.quad.html#scipy.integrate.quad\n", "\n", "val, err = quad(g, ga, gb)\n", "\n", "print(f\"Das Ergebnis unserer numerischen Berechnung mit der vorhandenen Funktion quad beträgt {val}.\")" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.12" }, "vscode": { "interpreter": { "hash": "cc412c9f3fa4474d685b397f63f04324611d10c22d923c0d092e91ec5ad6d5ee" } } }, "nbformat": 4, "nbformat_minor": 5 }