{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Runge Kutta módszerek (Explicit)\n", "\n", "\n", "Forrás: https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods\n", "\n", "### Kezdetiérték feladat (KDE):\n", "$\\dot{y}(t) = f(y,t)$ \n", "\n", "### Runge Kutta módszerek: \n", "Közelítsük a megoldást diszkrét időpontokban ($t_n = nh$): \n", "$y_n \\approx y(t_n)$ \n", "\n", "$\\begin{align}\n", "y_{n+1} = y_{n} + h \\sum_{i=1}^s b_i k_i,\n", "\\end{align}$ \n", "ahol \n", "$\\begin{align}\n", "k_1 &= \\left(t_n,y_n\\right), \\\\\n", "k_2 &= \\left(t_n+c_2 h,y_n+\\left(a_{21}k_1\\right) h\\right), \\\\\n", "k_3 &= \\left(t_n+c_3 h,y_n+\\left(a_{31}k_1+a_{32}k_2\\right)h\\right), \\\\\n", "&\\;\\vdots \\\\\n", "k_s &= \\left(t_n+c_s h,y_n+\\left(a_{s1}k_1+\\ldots+a_{s,s-1}k_{s-1}\\right)h\\right). \n", "\\end{align}$ \n", "\n", "Az $a_i$, $b_i$ és $c_i$ paramétereket az úgynevezett *Butcher tábla* (*Butcher tableau*) tartalmazza:\n", "\n", "$\\begin{equation}\n", "\\renewcommand\\arraystretch{1.2}\n", "\\begin{array}\n", "{c|cccc}\n", "0\\\\\n", "c_2 & a_{2,1}\\\\\n", "c_3 & a_{3,1} & a_{3,2} \\\\\n", "\\vdots & \\vdots& & \\ddots \\\\\n", "c_s & a_{s,1} & a_{s,2} & \\dots & a_{s,1}\\\\\n", "\\hline\n", "& b_1 & b_2 & \\dots & b_{s-1} & b_s\n", "\\end{array}\n", "\\end{equation}$\n", "\n", "#### Egyszerű esetek:\n", "Euler módszer: \n", "$\\begin{equation}\n", "\\renewcommand\\arraystretch{1.2}\n", "\\begin{array}\n", "{c|cccc}\n", "0\\\\\n", "\\hline\n", "& 1\n", "\\end{array}\n", "\\end{equation}$ \n", "$\\begin{align}\n", "y_{n+1} &= y_{n} + h \\cdot 1 \\cdot k_1,\\\\\n", "k_1 &= f(t_n, y_n).\n", "\\end{align}$\n", "\n", "Negyedrendű Runge-Kutta módszer: \n", "$\\begin{equation}\n", "\\renewcommand\\arraystretch{1.2}\n", "\\begin{array}\n", "{c|cccc}\n", "0\\\\\n", "\\frac{1}{2} & \\frac{1}{2}\\\\\n", "\\frac{1}{2} &0 &\\frac{1}{2} \\\\\n", "1& 0& 0& 1\\\\\n", "\\hline\n", "& \\frac{1}{6} &\\frac{1}{3} &\\frac{1}{3} &\\frac{1}{6} \n", "\\end{array}\n", "\\end{equation}$ \n", "$\\begin{align}\n", "y_{n+1} &= y_{n} + \\frac{h}{6} \\left(k_1+2k_2 + 2k_3 + k_4\\right),\\\\\n", "k_1 &= \\left(t_n,y_n \\right), \\\\\n", "k_2 &= \\left(t_n+\\frac{1}{2} h, y_n+\\frac{1}{2}k_1 h\\right), \\\\\n", "k_3 &= \\left(t_n+\\frac{1}{2} h, y_n+\\frac{1}{2}k_2 h\\right), \\\\\n", "k_4 &= \\left(t_n+ h, y_n+k_3 h\\right). \n", "\\end{align}$ \n", "\n", "#### Megjegyzés: Adaptív Runge-Kutta módszerek\n", "$\\begin{equation}\n", "\\renewcommand\\arraystretch{1.2}\n", "\\begin{array}\n", "{c|cccc}\n", "0\\\\\n", "c_2 & a_{2,1}\\\\\n", "c_3 & a_{3,1} & a_{3,2} \\\\\n", "\\vdots & \\vdots& & \\ddots \\\\\n", "c_s & a_{s,1} & a_{s,2} & \\dots & a_{s,1}\\\\\n", "\\hline\n", "& b_1 & b_2 & \\dots & b_{s-1} & b_s\\\\\n", "& b_1^\\ast & b_2^\\ast & \\dots & b_{s-1}^\\ast & b_s^\\ast\n", "\\end{array}\n", "\\end{equation}$ \n", "Hibafüggvény ami alapján eldöntjük, hogy kell-e felezni/duplázni $h$ értékét: \n", "$e_{n+1} = y_n+1-y^\\ast_n+1 = h \\sum_{i=s}^s\\left(b_i-b_i^\\ast\\right)k_i$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Implementáció" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "ExecuteTime": { "end_time": "2020-02-21T14:06:21.062475Z", "start_time": "2020-02-21T14:06:20.749695Z" } }, "outputs": [], "source": [ "# Számításokhoz és az eredmények megjelenítéséhez szükséges csomagok behívása\n", "import numpy as np\n", "import matplotlib.pyplot as mpl\n", "# Segédfüggvény\n", "def getzero(ll):\n", " return np.zeros(len(ll))\n", "\n", "# Runge Kutta lépés osztály\n", "class RKstep:\n", " def __init__(self,a,b,c,h):\n", " self.a = a\n", " self.b = b\n", " self.c = c\n", " self.h = h\n", " self.nk = len(c)\n", "\n", " def dy(self,f,t,y): # f(t_n,y_n) -> k\n", " ks=[getzero(y) for i in range(self.nk)]\n", " for (i,ci) in enumerate(self.c):\n", " ks[i]= f(t+ci*self.h,y + self.h*sum(self.a[i][j]*ks[j] for j in range(i)))\n", " return self.h*sum(k*self.b[i] for (i,k) in enumerate(ks))\n", "\n", "# RK4 módszer\n", "def RK4(h):\n", " a=[[],[1/2],[0,1/2],[0,0,1/2]]\n", " b=[1/6,1/3,1/3,1/6]\n", " c=[0,1/2,1/2,1]\n", " return RKstep(a,b,c,h)\n", "\n", "# Euler módszer\n", "def Euler(h):\n", " a=[[]]\n", " b=[1.]\n", " c=[0.]\n", " return RKstep(a,b,c,h)\n", "\n", "\n", "# Közönséges differenciálegyenlet kezdetiérték problémájának megadása/megoldásfüggvénye osztályként\n", "class ODEProblem:\n", " def __init__(self,f,y0,tspan):\n", " self.f = f\n", " self.y0 = y0 if type(y0)==np.ndarray else np.array(y0)\n", " self.tspan = tspan\n", " \n", " def solve(self,method):\n", " t0,T = self.tspan\n", " ts = np.arange(t0,T+100*np.finfo(float).eps,method.h)\n", " ys=[self.y0 for (i,t) in enumerate(ts)]\n", " for i in range(len(ys)-1):\n", " ys[i+1] = ys[i] + method.dy(self.f,ts[i],ys[i])\n", " return (ts,ys)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Alkalmazás\n", "\n", "Pl.: \n", "$\\dot{y}(t) = y(t), \\quad y(0) = y_0, \\quad t \\in \\left[0,T\\right]$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "ExecuteTime": { "end_time": "2020-02-21T14:06:21.429621Z", "start_time": "2020-02-21T14:06:21.064075Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxU5dXA8d9JSEjCDgk7MZFVQNaA+Lph1YBUQ61LiVXBpaiFam3pK2pft9YuahUVKsUFlGIQBWSQLWyKKCAJhJ2wBhICZAGykD3zvH9kwBAmZMhMMkvO9/PJZ2bufe7cc3P15PLc555HjDEopZTyXX7uDkAppVTd0kSvlFI+ThO9Ukr5OE30Sinl4zTRK6WUj2vk7gDsCQ0NNREREe4OQymlvEZiYmKWMSbM3jqPTPQREREkJCS4OwyllPIaInKkunXadaOUUj5OE71SSvk4TfRKKeXjPLKP3p7S0lLS0tIoKipydygeISgoiM6dOxMQEODuUJRSHs5rEn1aWhrNmjUjIiICEXF3OG5ljCE7O5u0tDQiIyPdHY5SysN5TddNUVERbdq0afBJHkBEaNOmjf7rRinlkBoTvYh0EZG1IrJHRHaJyNO25a1FZKWI7Le9tqpm+5EikiwiB0RksjPBapL/if4ulFKOcuSKvgz4ozHmKmAYMEFEegOTgdXGmO7AatvnC4iIPzANuB3oDcTatlVKKQWUlFl5ZfEujmSfrbN91JjojTHHjTFbbO/zgD1AJ2A08Imt2SfAL+xsPhQ4YIw5ZIwpAebatvNK/v7+DBgwgL59+3LvvfdSUFBQq++ZNWsW6enpLo5OKeWNvkhMZeb3KRzKdGOir0xEIoCBwCagnTHmOFT8MQDa2tmkE5Ba6XOabZm97x4vIgkikpCZmXk5YdWb4OBgkpKS2LlzJ4GBgUyfPr1W36OJXikFUFxWztQ1BxgY3pLhPe1WL3AJhxO9iDQF5gO/N8bkOrqZnWV2p7QyxswwxkQZY6LCwurugF3lhhtu4MCBAwC89dZb9O3bl759+zJlyhQAUlJSuOqqq/jNb35Dnz59iI6OprCwkC+//JKEhAR+/etfM2DAAAoLC1m9ejUDBw7k6quv5pFHHqG4uBiAyZMn07t3b/r168ekSZPcdqxKqbox98dUjucU8cfbetbpfTeHhleKSAAVSX6OMWaBbfFJEelgjDkuIh2ADDubpgFdKn3uDDh9KfvK4l3sTnf0b41jendszkt39nGobVlZGcuWLWPkyJEkJiYyc+ZMNm3ahDGGa665hptuuolWrVqxf/9+4uLi+OCDD7jvvvuYP38+DzzwAFOnTuXNN98kKiqKoqIixo0bx+rVq+nRowcPPfQQ77//Pg899BALFy5k7969iAhnzpxx6fEqpdyrqLScaWsPMDSyNdft/gHeiYfoaIiJcfm+HBl1I8BHwB5jzFuVVlmAsbb3Y4FFdjbfDHQXkUgRCQTG2LbzSoWFhQwYMICoqCjCw8N59NFHWb9+PXfddRdNmjShadOm/PKXv+S7774DIDIykgEDBgAwePBgUlJSLvrO5ORkIiMj6dGjBwBjx45l3bp1NG/enKCgIB577DEWLFhASEhIvR2nUqru/XfjETLyivljSCZyfyxMmwaxsWBxfYp05Ir+OuBBYIeIJNmWPQ/8A5gnIo8CR4F7AUSkI/ChMWaUMaZMRCYCKwB/4GNjzC5ng3b0ytvVzvXRV3apydUbN258/r2/vz+FhYUXtalu+0aNGvHjjz+yevVq5s6dy9SpU1mzZk0tI1dKeZKzxWW8/81Bru8WyjWbZsO5gR0FBRAf7/Kr+hoTvTFmPfb72gFusdM+HRhV6fNSYGltA/R0N954I+PGjWPy5MkYY1i4cCGzZ8++5DbNmjUjLy8PgF69epGSksKBAwfo1q0bs2fP5qabbiI/P5+CggJGjRrFsGHD6NatW30cjlKqHnyyIYXssyX8IboHtI2GmTMrknxISEX3jYt5TQkETzVo0CDGjRvH0KFDAXjssccYOHCg3W6ac8aNG8cTTzxBcHAwGzZsYObMmdx7772UlZUxZMgQnnjiCU6dOsXo0aMpKirCGMPbb79dT0eklKpLeUWlzFh3iJt7hjEovBWEx0BcXMWVfB310culuh7cJSoqylSdeGTPnj1cddVVborIM+nvRCnv886q/by9ah+WidfRr3NLl32viCQaY6LsrfOaWjdKKeXtcgpK+XD9IW7r3c6lSb4mmuiVUqqefLj+EHlFZfzhth71ul9N9EopVQ9OnS3h4/WH+fnVHbiqQ/N63bcmeqWUqgf/WXeQgtJyfn9r93rftyZ6pZSqY5l5xXz6wxFG9+9I93bN6n3/muiVUqqOvf/NQUrKrTx9a/32zZ+jif4yVC5TfOedd56vP5OSkkLfvn3Pt/vggw8YNGgQp0+fPr/szTffRETIysqq97iVUu5zIqeI/246wi8HdiIytIlbYtBEfxkqlylu3bo106ZNu6jN7Nmzee+994iPj6dVq4pJt1JTU1m5ciXh4eH1HbJSys2mrT2A1Wp46pb675s/RxN9LV177bUcO3bsgmXz5s3jH//4B/Hx8YSGhp5f/swzz/D666/r9H9KNTBppwuYu/ko9w3pQpfW7itM6NuJ3mKBiRNdXg2uvLyc1atXE1PpUeUjR44wceJE4uPjad++faUQLHTq1In+/fu7NAallOebtvYAgjDxZvfWqvLdRG+xVJT8dGHpz3Nlitu0acOpU6e47bbbzq8LCwsjPDycefPmnV9WUFDAa6+9xquvvur0vpVS3iUl6yxfJKQRO7QLHVsGuzUW30308fEXl/500rk++iNHjlBSUnJBH31ISAjLli1j+vTpzJkzB4CDBw9y+PBh+vfvT0REBGlpaQwaNIgTJ044HYtSyrP9Y9leAhv5McHNV/Pgy9Uro+uu9GeLFi149913GT16NE8++eT55WFhYSxfvpzhw4cTGhrKiBEjyMj4aeKtiIgIEhISLui/V0r5no2Hslm+6wSTonvQtnmQu8Px4Sv6GFvpzwkTKl5dXPpz4MCB9O/fn7lz516wPDIyEovFwiOPPMKmTZtcuk+llOezWg1/XbKbji2CeCxnT53cJ7xcWqbYi+nvRCnP82ViGpO+2MY73a2MfmrMT70KdXDBWZlTZYpF5GMRyRCRnZWWfS4iSbaflEpTDFbdNkVEdtjaJdhro5RSvqKgpIw3VuxlQJeWxGxd4fL7hLXlSNfNLGBk5QXGmF8ZYwYYYwYA84EFl9j+Zltbu39plFLKV0z/9hAnc4v5vzuuQqKjK67kocb7hJZkCxOXTsSSXDddPI7MGbtORCLsrZOKJ4DuA37m2rCUUsq7HM8pZMa6g9zRrwODr2gNVzg2RaAl2ULs/FgKSguYmTSTuLvjiOnp2i4eZ2/G3gCcNMbsr2a9AeJFJFFExl/qi0RkvIgkiEhCZmamk2EppVT9emN5MlYDz47s9dPCmBiYOvWSffPxB+MpKK3o4ikoLSD+oOu7eJxN9LFA3CXWX2eMGQTcDkwQkRura2iMmWGMiTLGRIWFhTkZllJK1Z9tqWdYsPUYj14fedmlDqK7RhMSULFNSEAI0V1dNxT8nFqPoxeRRsAvgcHVtTHGpNteM0RkITAUWFfbfSqllKcxpmI4ZWjTQH47vOtlbx/TM4a4u+OIPxhPdNdol3fbgHMPTN0K7DXGpNlbKSJNAD9jTJ7tfTTg1bUA/P39ufrqq89/HjNmDJMnT662/axZs0hISGDq1Kn1EZ5Syg2W7TzB5pTT/O2uq2kWFFCr74jpGVMnCf6cGhO9iMQBw4FQEUkDXjLGfASMoUq3jYh0BD40xowC2gELbRUbGwGfGWOWuzb8+nWuBEJdKSsro1Ej331YWSlfU1xWzt+X7aFX+2b8akgXd4dTLUdG3cRWs3ycnWXpwCjb+0NAgyjZWLm0QUJCApMmTeKbb765oE1mZiZPPPEER48eBWDKlClcd911vPzyy6Snp5OSkkJoaCifffaZG45AKVUbs75PIfVUIf999Br8/Ty3DLlPXz5aki0u7fc6V73ynOeee45f/epXDm379NNP88wzz3D99ddz9OhRRowYwZ49ewBITExk/fr1BAe7t8KdUspxWfnFTF1zgJ/1asv13T27fpXPJvq6GJvqTNfNqlWr2L179/nPubm55OXlARATE6NJXikvM2XVPgpLy3l+lOeXIfHZRG9vbGpd3exo1KgRVqsVgKKiIrttrFYrGzZssJvQmzRxzzySSqna2Xcyj882HeWhayPo1rapu8Opkc9Wr6yPsannREREkJiYCMD8+fPtxxMdfcHom7q8qauUqlt/XbKHpo0b8bQb54G9HD6b6M+NTZ0wZILLHik+10d/7ufc0MqXXnqJp59+mhtuuAF/f3+727777rskJCTQr18/evfuzfTp052ORylV/1bvOcm6fZk8dUt3WjUJdHc4DtEyxV5MfydK1a+CkjJue2sdTRr78/XvbiCwkedcKztVplgppVSFKav2c+xMIX87vJLApV+7OxyHaaJXSikH7ErP4aN1B4nduZKo916D2Fi3zxzlKK9K9J7YzeQu+rtQqv6UWw3PL9xJK2sJk1d9WLHQzZOJXA6vSfRBQUFkZ2drgqMiyWdnZxMU5P5Jh5VqCOZsOsK21DP8X88AWvjbclANk4lA3U8o4iivGUffuXNn0tLS0Fr1FYKCgujcubO7w1DK553MLeL15cnc0D2UmEeGQljNk4lA/Uwo4iivSfQBAQFERka6OwylVAPzyuJdlJZb+esv+iIiFcndgUm+6/OhzZp4TdeNUkrVt9V7TrJ0xwmeuqU7V7S5vCfY6/OhzZp4zRW9UkrVp4KSMl5ctIvubZvymxuuvOzt62NCEUdpoldKKTvesY2Z/+KJa2v9YFRdTyjiKO26UUqpKnan5/Lh+sPEDu3CkIjW7g7HaTUmehH5WEQyRGRnpWUvi8gxEUmy/YyqZtuRIpIsIgdEpPo595RSykNUjJnfQcvgAJ4d2cvd4biEI1f0s4CRdpa/bYwZYPtZWnWliPgD04Dbgd5ArIj0diZYpZSqa59tOkJS6hn+747etAzxjqJlNakx0Rtj1gGnavHdQ4EDxphDxpgSYC4wuhbfo5RS9eLcmPnru4UyekBHd4fjMs700U8Uke22rp1WdtZ3AlIrfU6zLbNLRMaLSIKIJOhDUUopd3h18W6KK4+Z9xG1TfTvA12BAcBx4F922tj7LVVbv8AYM8MYE2WMiQoLC6tlWEopVTsrd59kyY7jPPWzbkSE+tasb7VK9MaYk8aYcmOMFfiAim6aqtKALpU+dwbSa7M/pZSqS9n5xTy3YDtXdWjO+Bu7ujscl6tVoheRDpU+3gXstNNsM9BdRCJFJBAYA3hHTU+lVINhjOG5BTvIPVvMlINLvKrOvKMcGV4ZB2wAeopImog8CrwuIjtEZDtwM/CMrW1HEVkKYIwpAyYCK4A9wDxjzK46Og6llKqV+VuOEb/7JJO++y893/tnjXXmPaUi5eXwmqkElVLK1dJOFzByynf0Lsgg7l/j8DfWihUTJsDUqRe1r1yRMiQgxK0VKavSqQSVUqoKq9Uw6YttGGP415AW+Afb5ne4RJ15exUpvYHWulFKNUgff3+YjYdO8frd/egypAsE1VxnPrprNDOTZp6/ondnRcrLoV03SqkGZ9/JPO54bz03dg/lg4eiLmvMvCXZ4hEVKau6VNeNXtErpRqUkjIrz3yeRLPGjfj7L/td9oNRnlKR8nJooldKNSjvrdnPrvRcpj8wmLBmjd0dTr3Qm7FKqQZjy9HTTFt7gLsHdWZk3/buDqfeaKJXSjUIBSVl/OHzJDq0COalmIZVSFe7bpRSDcLfl+4lJbuAuN8Mo3lQgLvDqVd6Ra+U8nnf7stk9sYjPHp9JNd2bePucOqdJnqllE87U1DC/365je5tm/KnET3dHY5baNeNUspnWa2GP87bxqmzJXw0dghBAf7uDskt9IpeKeWz/rPuEKv3ZvBCwS76Jn7r7nDcRhO9UsonbTqUzZvL9/Dz/RsYO+V/fbIqpaM00SulfE5GXhET47YSbi3gH1+/VTHdXUFBRS0bO85VpZy2eRqx82N9LtlroldK+ZRyq+HpuCTyikp5f2AQzRrZShz4YFVKR+nNWKWUT3l75T42HMrmjXv60SuqCzTx3aqUjtLqlUopn7E2OYOHZ27mvqjOvH5P/8va1lOrUjrKqeqVIvIxcAeQYYzpa1v2BnAnUAIcBB42xpyxs20KkAeUA2XVBaGUUs46dqaQZz5Polf7Zrw6uu9lb++NVSkd5Ugf/SxgZJVlK4G+xph+wD7guUtsf7MxZoAmeaVUXSkpszJhzhbKyg3vPzC4wY6Xr06Nid4Ysw44VWVZvG3yb4CNQOc6iE0ppRzyt6V7SEo9wxv39CMytIm7w/E4rhh18wiwrJp1BogXkUQRGX+pLxGR8SKSICIJmZmZLghLKdUQLNl+nFk/pPDwdRHcfnUHd4fjkZxK9CLyAlAGzKmmyXXGmEHA7cAEEbmxuu8yxswwxkQZY6LCwsKcCUsp1UAcyszn2fnbGRjekuduv8rd4XisWid6ERlLxU3aX5tqhu4YY9JtrxnAQmBobfenlFKVFZaU89s5WwjwF6beP4jARvpYUHVq9ZsRkZHAs0CMMaagmjZNRKTZufdANLCztoEqpdQ5xhj+9OU2kk/m8favBtCpZbC7Q/JoNSZ6EYkDNgA9RSRNRB4FpgLNgJUikiQi021tO4rIUtum7YD1IrIN+BFYYoxZXidHoZRqUKas2s/X24/zbNFehidvvGRbX65h4yh9YEop5VUWJR3j6blJ3Lt7La8v/hcSEgJxcXafej1Xw+bcE69xd8f57Fj5Sz0wpZ1aSimvkXjkNH/6cjtDy0/x2pJ3aixW5us1bByliV4p5RXSThfw+OwEOrQIYvqwFgQGBVasuESxsuiu0YQEhFQ088EaNo7SomZKKY+XV1TKo7MSKC6zMnf8EFq3bQoBNRcri+kZQ9zdcV5dw8YVtI9eKeXRyq2Gxz7ZzLr9Wcx6eAg3dNfnbOzRPnqllNd6bcke1iZn8kpMH03ytaSJXinlseZsOsLH3x/m4esieGDYFe4Ox2tpoldKeaT1+7N4cdEuhvcM488/7+3ucLyaJnqllMc5kJHPk3MS6RrWhPdiB+LvJ+4OyatpoldKeZTTZ0t49JPNBPr78dHYITQLCnB3SF5PE71SymOcLS5j3NvxHM/KY0ZkIV1ah1TbVksbOE4TvVLKIxSVljP+reXszLUydcHfGPx4LFjsJ/FzpQ2mbZ5G7PxYTfY10ESvlHK7snIrT8Vt5fsc4fWl7xB9YJOWNnAhTfRKKbeyWg2TF+wgfvdJXoywcvdhWzVKLW3gMloCQSnlNsYY/rpkD18mpvH0Ld155LYe0FFLG7ialkBQSrnNu6v389bKfYz7nwheurM3IjqMsra0BIJSyuPM+v4wb63cxy8HdeLFOzTJ1yVN9EqperdwaxovL97NrVe14/W7++GnD0TVKUemEvxYRDJEZGelZa1FZKWI7Le9tqpm25EikiwiB0RksisDV0p5p5W7TzLpi+1ce2Ubpt4/kEb+er1Z1xz5Dc8CRlZZNhlYbYzpDqy2fb6AiPgD04Dbgd5ArIhowQqlGrANB7OZ8NkW+nZszgdjowgK8Hd3SA1CjYneGLMOOFVl8WjgE9v7T4Bf2Nl0KHDAGHPIGFMCzLVtp5RqgJJSz/CbTxO4onUIMx8eStPG9gf96ROvrlfbfzO1M8YcB7C9trXTphOQWulzmm2ZXSIyXkQSRCQhMzOzlmEppTxR4pFTPDh9Pa0Kc5ndJYfWTQLtttMnXutGXXaO2bu7Uu1YTmPMDGNMlDEmKixMJxdQyldsPJTNg//5gdCs43w+/be0H1d9aQN94rVu1DbRnxSRDgC21ww7bdKALpU+dwbSa7k/pZQXWr8/i3Ezf6Rj2Vk+n/MsHfOyLlnaQJ94rRu1TfQWYKzt/VhgkZ02m4HuIhIpIoHAGNt2SqkGYO3eDB75ZDMRbZowd2gwbU1xxYpLlDY498TrhCETiLs7Tp94dZEaSyCISBwwHAgVkTTgJeAfwDwReRQ4Ctxra9sR+NAYM8oYUyYiE4EVgD/wsTFmV90chlLKk8TvOsGEz7bQs30zZj9yDa2aBEJgzaUNoCLZa4J3LS2BoJRyqSXbj/P03K307dSCTx4ZSotgnTikPmgJBKVUvVi4NY3fxW1hYHhLZj+qSd5TaPVKpZRLzNucyrMLtjMssg0fjYsiJFDTi6fQM6GUctrsjUf4v692cmOPMGY8OFifePUwmuiVUrVmjGHG9K/5+xE/bmllmHaJJG9Jtmj9eDfRPnqlVK2UWw0vvVOR5H++Zx3v//V+gpYtsdtWn3h1L030SqnLVlhSzhP/TeTTE36M3zSf9yxvEJifp3O8eihN9Eqpy5KdX0zsBxtZteckr0Raef7Hz/HD6ByvHkz76JVSDjucdZZxM3/kRE4R0x8YzIg+7aGDzvHq6fSBKaWUQ7YcPc1jn1T8f/nh2CgGhdudb0i5yaUemNIreqVUjVbsOsFTcVtp3yKITx4eSkRoE3eHpC6DJnql1CXN+v4wr3y9m/6dW/LR2CjaNG3s7pDUZdJEr5Syy2o1/H3ZHj747jC39W7Hu2MGEhx48Rh5HR/v+TTRK6UukldUyqR3V7DilDC2vZUXHxiMv9/FcwmdGx9fUFrAzKSZWlrYQ+nwSqXUBQ5m5vOLf65gVZaVF1fN4OUXxuD/9WK7bXV8vHfQRK+UOi9+1wlGT/2eM2dL+O/cP/NIogXRGaG8nnbdKKWwWg1TVu3j3TUH6Ne5BdPbFdNx2sGKlQ7MCKV99J6t1oleRHoCn1dadCXwojFmSqU2w6mYZvCwbdECY8yrtd2nUsr1cgpLeebzJNbszeDewZ35yy/6VhQma6wzQvmKWid6Y0wyMABARPyBY8BCO02/M8bcUdv9KKXqzr6TeTw+O5HUUwX85Rd9eeCacERsN11jYi6Z4JX3cFXXzS3AQWPMERd9n1Kqji3dcZxJX2yjSeNGxI0fxpCI1u4OSdURV92MHQPEVbPuWhHZJiLLRKSPi/anlKqlcqvhn8v38ts5W+jVvhlf/+56u0nekmxh4tKJWlLYBzhd60ZEAoF0oI8x5mSVdc0BqzEmX0RGAe8YY7pX8z3jgfEA4eHhg48c0X8cKOVqJ3KKeGbaKjbkCve3s/LS70bRuJH9h6DOjY8PCQjR8fFeoK4nB78d2FI1yQMYY3KNMfm290uBABEJtfclxpgZxpgoY0xUWFiYC8JSSlW2YtcJRr6xim1Zxby+dAp/+/MYGi+1P1GIjo/3La5I9LFU020jIu3FdmdHRIba9pftgn0qpRxUWFLO8wt38PjsRLoU5/L1rKe4b8cq0PHxDYZTN2NFJAS4DXi80rInAIwx04F7gCdFpAwoBMYYT6yLrJSP2p2ey1Nzt3IgI5/Hb7ySP5bsJ/DfZypW6vj4BkPr0Svlg4wxzPw+hX8s20vLkADeum8A13e39ZpaLA6Nj1feRevRK9WAZOUXM+mLbXyTnMmtV7Xln3f3u6C0sKUnxDeC6K6gab5h0ESvlA/5JjmDSV9sJ7eolFdH9+HBYVf89AAUWm2yodJEr5QPOFtcxhvvL2PWCT96hhjmTLyJnu2bXdTO3mgaTfS+T6tXKuXlvtufSfTflvPJcRiXYGHRmw/Q88e1dtvqaJqGSa/olfJSOQWl/HXJbr5ITONK61nmffYqQ47trlgZH2/3RquOpmmYNNEr5YVW7DrBn7/ayamzJfx2eFeeKtpH0LSUipWXGDYJWm2yIdJEr5QXycwr5mXLLpbsOE7vDs2ZOW4IfTu1AHpBnGNlhVXDo4leKS9gjOGrpGO8sng3BcXl/GlET8bfeCUB/j/dZtNhk6o6muiV8nDHvrDwwqZsvmnUlkHhLXn9nn50a3vhiBodNqkuRUfdKOWhikrLmfrvxdy6sZRN1ua8uG4mX3TKvijJgxYhU5emiV4pD2OMYeXuk0S/vY43j/ox/FAi8R/9lkc2zMd/pRYhU5dPu26U8iCHMvN5ZfFuvt2XSbe2Tflv73Kuf++dikqTWoRM1ZImeqU8QH5xGe+t2c/H6w8T1MifP//8Ksb+TwQB/n5Y5PfE77QQ3TeGGJ2kW9WCJnql3MgYw6KkdP62dA8ZecXcM7gzz47sRViziiJklmQLsalTKAgpYGbqIeKSr9Fkri6bJnql3GTnZxZe3ppDgn9r+nduwX8eHMzA8FYXtNHaNMoV9GasUvUsJessE99YzB3b/Tlc5Mc/V7/Pwi6nLkryoDdZlWvoFb1S9SQjt4h3Vu/n882pBJRbmbhhHuM3LaB5SQGsjIDRWptG1Q1npxJMAfKAcqCs6uwmtvli3wFGAQXAOGPMFmf2qZS3ySks5T/fHuTj7w9TVm64/5pwJpYeou3UBVBy6dE0oDdZlfNccUV/szEmq5p1twPdbT/XAO/bXpXyeUWl5XzyQwr//uYgOYWljB7QkT/c1oMr2jQB+mIpdGw0jVLOquuum9HAp7YJwTeKSEsR6WCMOV7H+1XKbcoWWfhy7W6mNO/LiRJheM8w/jSiJ306tjjfRkfTqPrk7M1YA8SLSKKIjLezvhOQWulzmm3ZRURkvIgkiEhCZmamk2EpVf9Kyqx8/sFibl2WweSgq+lwZD9z+5Qz6+GhFyR50JIFqn45m+ivM8YMoqKLZoKI3FhlvdjZxtj7ImPMDGNMlDEmKiwszMmwlKo/57pohr+xlmcP+tG06Cwz5v+FBZ/+gWGbtGSBcj+num6MMem21wwRWQgMBdZVapIGdKn0uTOQ7sw+lfIU+cVlzNl4hA++O0xWfjFRV7TitY5nyV32v1i6l1AUEkiMlixQHqDWiV5EmgB+xpg82/to4NUqzSzARBGZS8VN2Bztn1feLqeglFk/pDDzh8OcKSjl+m6hTPzZQK6JbM3ifYu5/26hwMBMEeJ6Vl8bXkfTqPrizBV9O2BhxQhKGgGfGWOWi8gTAMaY6akjKEYAAA5TSURBVMBSKoZWHqBieOXDzoWrlPtkfWnhox9SmB18Jfnlwq1XtWXCzd0ueNAp/mA8BaYYgAJTrE+yKo9Q60RvjDkE9LezfHql9waYUNt9KOUJ9hzPZWbct3yVXk5pQCSj9m5gwj3X0Dt2yEVto7tGMzNpJgWlBdr3rjyGPhmrlB3lVsOavRl8vP4wGw5lE2ys3LtjFY8kLKLrqWPQYQLE6pOsyjtooleqkvziMr5ISGXWDykcyS6gY4sgJt/eizFZu1i3/EPevqaE6GPV32QF7XtXnkcTvVIWC0dXfsesbjfyxalA8orLGHxFK/53RC9G9GlHI38/LMm7iHXwJqtSnkYTvWqwysqtfDv7a+KWbmF15I34p1m5o10JDz92A/27tLygrd5kVd5ME71qcNJOFzBvcyrzEtI4ketPaNtuTNgwjwe3LqXd2FjocudF2+hNVuXNNNGrBqG03Mqq3SeJ25zKd/srSmzc1COMlzucpfCd8azuWMimK/UBJ+WbNNEr32WxcHjleuZ2v4H5uUFk5ZfQoUUQT/2sO/cN6UKnlsFYki08kGrVB5yUT9NEr3xOTmEpy/67nIXf7mVT55vwTyvnZ6HFxI4bwk092uLv91MJJu17Vw2BJnrlE4rLylm7N4OvtqazZm8GJeV+tG3Rgp4nPuWBpFU8eH0s9Lrjou207101BJroldeyWg0/ppxiUdIxlmw/Tm5RGaFNG/PrYeG0Ob6IZ4+8QEFbWN8ZWvRsYbdLRvveVUOgiV55FWMMe+cuZtHmFCxNryS9RAgJ9Gdkn/aMHtiJ67q2oZG/HxOXvk+BrXxeQQDEt8nRvnfVYGmiVx7PGMO2tByW7zzB8k0HSCnyRwIi6J2WyLPXDeK2B0cQEnjhf8raJaPUTzTRK49UbjUkHjnNsp3HWbHzBOk5RTTyE7qVpZLnv4gzQRvI6JDL47tGERJ48bh37ZJR6iea6JXHKP3Kwsa1iSyLiCL+bDBZ+cUENvLjxu5h/CG6J7de1ZY/z7yLFVkrACgIhPiOOhxSqZpooldulZ1fzLf7MlmzeiurTxRSGDyUwGOF3NYORt45mJt7taVp45/+M42+9XFmfr6aAlNMiDQm+tbH3Ri9Ut5BE72qV8YYdqXnsmZvBmv2ZrAt7QzGQDNTxKnADeQFbET8tjLi+K3c2X/JRdvH9Iwh7lfztEtGqcugiV7VufyFFtZ/s421VwxgbWEwGXnFiEC/zi35/S09uLlXGP/56tf8O3tpxQb+EN9Ju2SUchVn5oztAnwKtAeswAxjzDtV2gwHFgGHbYsWGGOqziurfExpuZVtqWdYfyCL7zftIzHHYA0eRNCxs9zSHm6+dyDDe4YR2rTx+W1G3PY4s7RLRqk64cwVfRnwR2PMFhFpBiSKyEpjzO4q7b4zxlz8SKLyGcYY9p3MZ/2BLH44kMXGQ9mcLSlHBLqUZ5Lv/y25AVvw99vNiOMjiBmsXTJK1Sdn5ow9Dhy3vc8TkT1AJ6Bqolc+xhjDkewCZnz2Id8eyyMvoA+5ZYEARIY24a5Bnbi+WyjDrmzDn2fexXdZ2iWjlDu5pI9eRCKAgcAmO6uvFZFtQDowyRizq5rvGA+MBwgPD3dFWMpFrFbD3hN5bE45xY8pp9h8+BQZecVAN8rlNKXlP/BIx85MeGg8nVoGX7CtjpJRyv2cTvQi0hSYD/zeGJNbZfUW4ApjTL6IjAK+Arrb+x5jzAxgBkBUVJRxNi5VeyVfWdix5kd+7HUNm0Pak5ByityiMgA6tAji2q5tOLR1GiuCVlIqR0Eg91hfOrV8+qLv0i4ZpdzPqUQvIgFUJPk5xpgFVddXTvzGmKUi8m8RCTXGZDmzX+U6xhjSc4rYevQ0SUfPsHXbIbafKaM05Fo4ClcGn+Tn/cIZEtGaIRGt6dwqGBHBUtqONclHKQ2AkFKI7lt9AtcuGaXcy5lRNwJ8BOwxxrxVTZv2wEljjBGRoYAfkF3bfSonWSzkx69m+zW3kNSlN1uPniEp9QyZeRX12AMb+dGp5AQ5jTZxNmA3Urabd3P+h5hf2rl5+uBrxM2G+J0WovvGEPPga/V9NEopBzlzRX8d8CCwQ0SSbMueB8IBjDHTgXuAJ0WkDCgExhhjtFumnhSUlLE7PZedx3LYuXEnPxxMJb3ZLbDLH3YlExnahOu7hTIwvCUDurSkV/vm/OE/d/KNozdPH3yNGDTBK+XpnBl1sx6QGtpMBabWdh/KMZbZL7Bk5wp6RN5Liz73ViT29FwOZuZz7s9qM1NEdvMsChr9gFiT+XdROLGTLupt05unSvkgfTLWy5RbDUdPFbD3eC57T+Sx9seNbM2LwF9eYUUKkLKHDi2C6NOxBXf060Dfji24unMLXv307p+ePAW+79SOWDvfrzdPlfI9mug92JyZL7N4bxIdOtxKQOcb2Hsij30n8ygqtQLgJxBstVLkf5gSv5WU+B3kgeJAPnhuw0XfdTlPnurNU6V8iyZ6d7BYID4eoqOx3nEn6TmFHMjI50BGPgczK153p2ZwtnwIyBA4AU0zD9M/sj33D72CXh2a0at9M7q3bcbKz18kNvnvFNhGv9zZ93m7u9QrdaUaLvHEe6NRUVEmISHB3WG4VEFJGYezzpKy/FuWLZ7L1rYd8PMLJ7tVVwqsP93qaBUSQLe2TTl5eDlJgdsplRRK/I7wZGFnpv5zh93vtsx+QUe/KNXAiUiiMSbK3jq9oncRy+wXWLZzKQO63U3E/zxWkdSzzla8Zp/lZG6xraU/9Pw1ZZKJ1aRyW+kObrv7frq1bUr3tk1pYyv0ZZm9hNjkRRQ1OjdO/bfV7ltHvyilLkUTfQ0qXy3f8eu/kplfTOqpAo7aflJPFZK0dxf7CnrgJ8NYdtAPDiYC0LpJoG0IYxiRoSFEhDbh89XPMSfnC4wUAtAsdBQPDHvuov3qOHWllKs02K4be90d5VZDRl4R6WcKSTtdyKo1i5iXuRMjbQm0tiNYOlFq/C/4nvbNgyjJ2cPRRocokxOUSjp3lQQz7WULLYIDLt5vsoXYz+87f1M07lfztL9cKeW0BtV1U11/tTGGU2dLOJ5ThGXRB/wrbSvWgAF8ubuId179jOLG7TmRU0SZtfIfvisIpDXlZFDil0aP8hQe/sWTdG4dQnjrEDq1DCYowB/L7DXEJr99/obofX2ft5vkQW+KKqXqn09d0Vtmv8ADez8Fv0EElYVyfbMb8G/dgxO5RZzMKaak3HpBe0M55ZJNW+tZbhx4HR1bBtOpVTAdWwbTuWUwW5e9zsMH/no+gcf1fL7aLhS9IaqUcqcGc0Ufv9OCNbAHoaXjsVLM3pw8rm4jDApvRfsWQbRvHkSHFkHs+34WLxz9F/kBpwkps/Lv3s8TM2bgRd/XfdxfCJptdSiB6w1RpZSn8rkr+ti9b1EUEEBQaR5xvfQKXCnVMFzqit6nEj1oAldKNUwNKtErpVRDdKlE71ffwSillKpfmuiVUsrHaaJXSikf51SiF5GRIpIsIgdEZLKd9SIi79rWbxeRQc7sTyml1OWrdaIXEX9gGnA70BuIFZHeVZrdDnS3/YwH3q/t/pRSStWOM1f0Q4EDxphDxpgSYC4wukqb0cCnpsJGoKWIdHBin0oppS6TM4m+E5Ba6XOabdnltlFKKVWHnEn09iYGrzoo35E2FQ1FxotIgogkZGZmOhGWUkqpypypdZMGdKn0uTOQXos2ABhjZgAzAEQkU0SOOBFbKJDlxPaexFeORY/Ds+hxeB5nj+WK6lY4k+g3A91FJBI4BowB7q/SxgJMFJG5wDVAjjHmeE1fbIwJcyIuRCShuifEvI2vHIseh2fR4/A8dXkstU70xpgyEZkIrAD8gY+NMbtE5Anb+unAUmAUcAAoAB52PmSllFKXw6kyxcaYpVQk88rLpld6b4AJzuxDKaWUc3z1ydgZ7g7AhXzlWPQ4PIseh+eps2PxyOqVSimlXMdXr+iVUkrZaKJXSikf59WJ3leKqjlwHMNFJEdEkmw/L7ojzpqIyMcikiEiO6tZ7y3no6bj8Jbz0UVE1orIHhHZJSJP22nj8efEwePw+HMiIkEi8qOIbLMdxyt22tTN+TDGeOUPFUM6DwJXAoHANqB3lTajgGVUPKE7DNjk7rhreRzDga/dHasDx3IjMAjYWc16jz8fDh6Ht5yPDsAg2/tmwD4v/X/EkePw+HNi+x03tb0PADYBw+rjfHjzFb2vFFVz5Di8gjFmHXDqEk284Xw4chxewRhz3BizxfY+D9jDxbWmPP6cOHgcHs/2O863fQyw/VQdDVMn58ObE72vFFVzNMZrbf/kWyYifeonNJfzhvPhKK86HyISAQyk4iqyMq86J5c4DvCCcyIi/iKSBGQAK40x9XI+nHpgys1cWlTNjRyJcQtwhTEmX0RGAV9RUePf23jD+XCEV50PEWkKzAd+b4zJrbraziYeeU5qOA6vOCfGmHJggIi0BBaKSF9jTOV7QXVyPrz5it6lRdXcqMYYjTG55/7JZyqeRg4QkdD6C9FlvOF81MibzoeIBFCRHOcYYxbYaeIV56Sm4/CmcwJgjDkDfAOMrLKqTs6HNyf680XVRCSQiqJqliptLMBDtjvZw3CwqFo9q/E4RKS9iIjt/VAqzlt2vUfqPG84HzXylvNhi/EjYI8x5q1qmnn8OXHkOLzhnIhImO1KHhEJBm4F9lZpVifnw2u7boyPFFVz8DjuAZ4UkTKgEBhjbLfoPYmIxFEx+iFURNKAl6i44eQ15wMcOg6vOB/AdcCDwA5bvzDA80A4eNU5ceQ4vOGcdAA+kYppWP2AecaYr+sjZ2kJBKWU8nHe3HWjlFLKAZrolVLKx2miV0opH6eJXimlfJwmeqWU8nGa6JVSysdpoldKKR/3/+ZY/xstKWN8AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Függvény és ODEProblem megadása\n", "foo = lambda t,y: y#*np.sin(4*t)\n", "y0=1.\n", "T=3.\n", "prob=ODEProblem(foo,[y0],(0.,T));\n", "\n", "# Lépésköz + numerikus megoldás\n", "h=0.1\n", "# Runge Kutta megoldás\n", "ts,ys=prob.solve(RK4(h))\n", "#Euler módszer megoldás\n", "tsE,ysE=prob.solve(Euler(h))\n", "\n", "# Ábrázolás\n", "mpl.plot(ts,np.exp(ts))\n", "mpl.scatter(ts,[y[0] for y in ys],color='r',s=10)\n", "mpl.scatter(tsE,[y[0] for y in ysE],color='g',s=10)\n", "mpl.legend(['Pontos','RK4','Euler'])\n", "mpl.show()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.7.3" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }