{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2 · Galerkin-Verfahren und Trapezregel\n",
    "\n",
    "Wir zerlegen das Intervall in gleich große Elemente und verwenden stückweise lineare Hutfunktionen $\\varphi_i$. Der Ansatz lautet\n",
    "\n",
    "$$u_h(x,t)=\\sum_i c_i(t)\\varphi_i(x).$$\n",
    "\n",
    "Nach Multiplikation mit Testfunktionen und partieller Integration entsteht\n",
    "\n",
    "$$M\\dot c(t)+\\kappa Kc(t)=0.$$\n",
    "\n",
    "$M$ ist die Massenmatrix und $K$ die Steifigkeitsmatrix. Für den Zeitschritt verwenden wir Crank–Nicolson, also die Trapezregel:\n",
    "\n",
    "$$(M+\\tfrac{\\Delta t\\kappa}{2}K)c^{n+1}=(M-\\tfrac{\\Delta t\\kappa}{2}K)c^n.$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys\n",
    "from pathlib import Path\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "project = Path.cwd()\n",
    "if not (project / 'galerkin_heat.py').exists() and (project.parent / 'galerkin_heat.py').exists():\n",
    "    project = project.parent\n",
    "if not (project / 'galerkin_heat.py').exists():\n",
    "    from urllib.request import urlretrieve\n",
    "    source = 'https://raw.githubusercontent.com/karkessler/dhbw-mathe3/main/galerkin/galerkin_heat.py'\n",
    "    urlretrieve(source, project / 'galerkin_heat.py')\n",
    "sys.path.insert(0, str(project))\n",
    "\n",
    "from galerkin_heat import exact_solution, l2_error, solve_heat_equation\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compare(elements=30, dt=0.004, final_time=0.1):\n",
    "    solution = solve_heat_equation(elements, dt, final_time)\n",
    "    reference = exact_solution(solution.x, final_time)\n",
    "\n",
    "    plt.figure(figsize=(8, 4))\n",
    "    plt.plot(solution.x, reference, '--', lw=2, label='exakt')\n",
    "    plt.plot(solution.x, solution.values[-1], 'o-', ms=3, label='Galerkin')\n",
    "    plt.xlabel('Ort x')\n",
    "    plt.ylabel('u(x,t)')\n",
    "    plt.title(f'{elements} Elemente, Δt={dt:g}, Fehler={l2_error(solution):.2e}')\n",
    "    plt.grid(alpha=0.3)\n",
    "    plt.legend()\n",
    "    plt.show()\n",
    "\n",
    "compare()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gitter und Zeitschritt verändern\n",
    "\n",
    "Mehr Elemente verbessern die räumliche Auflösung. Ein kleineres $\\Delta t$ reduziert den Zeitdiskretisierungsfehler."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    from ipywidgets import FloatLogSlider, IntSlider, interact\n",
    "    interact(\n",
    "        compare,\n",
    "        elements=IntSlider(value=30, min=6, max=80, step=2),\n",
    "        dt=FloatLogSlider(value=0.004, base=10, min=-3.0, max=-1.3, step=0.1),\n",
    "        final_time=FloatLogSlider(value=0.1, base=10, min=-2.0, max=-0.5, step=0.1),\n",
    "    )\n",
    "except ImportError:\n",
    "    print('Für Regler bitte ipywidgets installieren.')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Einordnung\n",
    "\n",
    "Das räumliche Verfahren ist ein Galerkin-Verfahren, verwendet hier aber bewusst einfache lokale Hutfunktionen. Dadurch bleibt das Beispiel übersichtlich und die Implementierung unabhängig vom historischen Programmlisting."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"},
  "language_info": {"name": "python", "version": "3"}
 },
 "nbformat": 4,
 "nbformat_minor": 5
}