P02b-ClassicalConstForceMotion.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Classical 1-D acceleration with The Heun Method\n",
    "-----------------------------------------------\n",
    "\n",
    "Suppose we have a constant force $F_0$ acting on a fixed mass for a period of time. At low speeds the momentum principle gives us:\n",
    "\n",
    "$$\\frac{dv}{dt} = \\frac{F_0}{m}$$\n",
    "\n",
    "and we also know that\n",
    "\n",
    "$$\\frac{dx}{dt} = v$$\n",
    "\n",
    "which form a set of two coupled differential equations to be solved for the position and velocity of the mass."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x105d7c5c0>]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAHDJJREFUeJzt3X28XFV97/HP1xBA5cmSUw0h4QiiFCyPRyC12hRKi9FLbi/wanxAodhcVFp4iV6Ba1HBitYX1AIq5pZcCLWABS8NFGsRAoESIieRpyRQAwVJDHIIkgR50MDv/rHXGYbJzJk958yex+/79ZpXZu+9Zu/fmg3nN3utvdZWRGBmZgbwunYHYGZmncNJwczMSpwUzMysxEnBzMxKnBTMzKzEScHMzEqcFKwlJF0q6a879fiSvijpH1sZUyeT9GZJSyRtlnRBu+Ox1nFSsNeQ9JikX0uaUrH+J5JC0mCOfZwo6c7ydRFxSkSc19xo8ys/vqRZkta2K5ZGNTveaueninnA08BOEXFGs45tnc9Jwar5L+CDowuSfhd4Q/vCsTbYA1gV4xjdKmmbAuKxFnFSsGquBD5atvwxYGF5AUk7S1ooaUTS45I+L+l1kn4HuBSYKek5Sc+m8pdL+nLZ5/9C0hpJz0haJGm3sm0h6RRJP5X0rKRvSlJlkJK2l/TC6FWNpP8taYukndLyeZK+UX58SW8EfgDsluJ7ruzY26Y6bZa0UtJQrS9I0n6Sbk7x/0LS2Wn9dpK+Ienn6fUNSdulbbMkrZV0hqSnJK2XdFLZPmdLWpWOv07SZ2rFK+lQSUvT97Ne0iWStq33HdY6PxV1uzyd8/+VyvxRznp9TtKTwP9N6+dIulfSJkmPSDq67L+dy1Lc69J5mVTru7YWiwi//Cq9gMeAPwIeBn4HmASsJfvlGMBgKrcQ+BdgR2AQ+E/g5LTtRODOiv1eDnw5vT+CrGniYGA74GJgSVnZAG4EdgFmACPA0TXiXQIcm97/O/AI8L6ybX9a5fizgLUV+/ki8CIwO9X5fODuGsfcEVgPnAFsn5YPS9vOBe4GfhsYAO4Czis77pZUZnI61vPAm9L29cB70vs3AQePEe8hwOHANun7Xw2cnuc7rHZ+qtSx9H01UK+vpfP5euBQYCNwFNmPz2nAPqn8/wO+A7wx7e/HwP9s93/7fmUvXylYLaNXC0eR/cFZN7oh/aqbC5wVEZsj4jHgAuCEnPv+MLAgIlZExEvAWWS/XAfLynw1Ip6NiJ8Bi4EDa+zrduAPUpPF/sBFaXl74F1kiSGvOyPipoh4maz+B9Qo9wHgyYi4ICJeTN/BsrK6nRsRT0XECPAlXvu9/CZt/01E3AQ8B7yjbNu+knaKiF9GxIpagUbE8oi4OyK2pO//O8AfVBTL+x3mUa9erwBfiIiXIuIF4GSyc3xzRLwSEesi4iFJbyZLhqdHxK8i4ing78j+e7IO4KRgtVwJfIjsV+XCim1TyH7pPl627nGyX4N57Fb+2Yh4DthQ8fkny94/D+xQY1+3k/1SPRh4ALiZ7I/j4cCaiNiQM6Zqx9y+Rvv4dLIrkmpeU7f0frey5Q0RsaXiOKN1O5bsD+bjkm6XNLNWoJLeLulGSU9K2gR8hey8jFWfWt9hHvXqNRIRL5Yt1/qO9iD7b2d9atZ6liyh/fYEYrMmclKwqiLicbIO59nA9ys2P032q3aPsnUzePVqol7n5M/LP5vazXct+3wj7iL7pf2nwO0RsSrFMpssYVQz0amBnwD2rLHtNXVLsfw8z04j4p6ImEP2B/J64Hujm6oU/zbwELB3ROwEnA1s1e9S61A5y5WrV6/KfT4B7FVlP08ALwFTImKX9NopIvYbR0xWACcFG8vJwBER8avylal55XvA30jaUdIewKeB0fv8fwHsXt7xWeEq4CRJB6bOyq8Ay1IzSEMi4nlgOfApXk0CdwGnUDsp/ALYVdLOjR4vuRGYKun01AG7o6TD0rargM9LGkgd4Ofw6vdSk6RtJX1Y0s4R8RtgE1mTTK14d0xlnpO0D/CJBuKvd36qabRel5Gd4yOV3YAwTdI+EbGerO/nAkk7pW17Saps+rI2cVKwmiLikYgYrrH5L4FfAY8CdwL/BCxI224FVgJPSnq6yn5/BPw1cB1Z5+peTKxN+XayJokfly3vSI3+hIh4iOyP3KOpCWO3auVqiYjNZH0t/42sieanwB+mzV8GhoH7yZqzVqR1eZwAPJaag04ha8evFe9nyJr3NgP/B7imgSqMeX5qaKheEfFj4CSy/oKNZOdk9Erjo8C2wCrgl8C1wNQG4rcCKcIP2TEzs4yvFMzMrMRJwczMSpwUzMysxEnBzMxKum7iqilTpsTg4GC7wzAz6yrLly9/OiIG6pXruqQwODjI8HCtuyTNzKwaSY/XL+XmIzMzK+OkYGZmJU4KZmZW4qRgZmYlTgpmZlbipGBmZiVOCmZm3WDpUjj//OzfAnXdOAUzs76zdCkceST8+tew7bZwyy0ws+aD+SbEVwpmZp3uttuyhPDyy9m/t91W2KGcFMzMOt2sWdkVwqRJ2b+zZhV2KDcfmZl1upkzsyaj227LEkJBTUfgpGBm1h1mziw0GYxy85GZmZUUlhQkbS/px5Luk7RS0peqlNlO0jWS1khaJmmwqHjMzKy+Iq8UXgKOiIgDgAOBoyUdXlHmZOCXEfE24O+ArxUYj5mZ1VFYUojMc2lxcnpFRbE5wBXp/bXAkZJUVExmZja2QvsUJE2SdC/wFHBzRCyrKDINeAIgIrYAG4Fdq+xnnqRhScMjIyNFhmxm1tcKTQoR8XJEHAjsDhwq6Z3j3M/8iBiKiKGBgbpPkzMzs3Fqyd1HEfEssBg4umLTOmA6gKRtgJ2BDa2IyczMtlbk3UcDknZJ718PHAU8VFFsEfCx9P444NaIqOx3MDOzFily8NpU4ApJk8iSz/ci4kZJ5wLDEbEIuAy4UtIa4BlgboHxmJlZHYUlhYi4Hzioyvpzyt6/CBxfVAxmZtYYj2g2M2uHFj0foVGe+8jMrNVa+HyERvlKwcys1Vr4fIRGOSmYmbVaC5+P0Cg3H5mZtVoLn4/QKCcFM7N2aNHzERrl5iMzMytxUjAzsxInBTMzK3FSMDOzEicFMzMrcVIwM7MSJwUzMytxUjAzsxInBTMzK3FSMDObqA6dBns8PM2FmdlEdPA02OPhKwUzs4no4Gmwx8NJwcxsIjp4GuzxcPORmdlEdPA02OPhpGBmNlEdOg32eLj5yMzMSpwUzMyspLCkIGm6pMWSVklaKem0KmVmSdoo6d70OqeoeMzMrL4i+xS2AGdExApJOwLLJd0cEasqyt0RER8oMA4zM8upsCuFiFgfESvS+83AamBaUcczM7OJa0mfgqRB4CBgWZXNMyXdJ+kHkvar8fl5koYlDY+MjBQYqZlZfys8KUjaAbgOOD0iNlVsXgHsEREHABcD11fbR0TMj4ihiBgaGBgoNmAzsz5WaFKQNJksIXw3Ir5fuT0iNkXEc+n9TcBkSVOKjMnMzGor8u4jAZcBqyPiwhpl3pLKIenQFM+GomIyM8ulh2Y9bVSRdx+9GzgBeEDSvWnd2cAMgIi4FDgO+ISkLcALwNyIiAJjMjMbW4/NetqowpJCRNwJqE6ZS4BLiorBzKxh1WY97aOk4BHNZmblemzW00Z5Qjwzs3I9Nutpo5wUzMwq9dCsp41y85GZmZU4KZiZWYmTgpmZlTgpmJlZiZOCmZmVOCmYmVmJk4KZmZU4KZhZ7+vjCe4a5cFrZtbb+nyCu0b5SsHMelu1Ce6sJicFM+ttfT7BXaPcfGRmva3PJ7hrlJOCmfW+Pp7grlFuPjIzsxInBTMzK3FSMDOzEicFMzMrcVIwM7MSJwUzMyspLClImi5psaRVklZKOq1KGUm6SNIaSfdLOrioeMzMrL4irxS2AGdExL7A4cCnJO1bUeZ9wN7pNQ/4doHxmFmv8AR3hSls8FpErAfWp/ebJa0GpgGryorNARZGRAB3S9pF0tT0WTOzrXmCu0LlulKQ9CZJ+0naU1LDVxeSBoGDgGUVm6YBT5Qtr03rKj8/T9KwpOGRkZFGD29mvcQT3BWq5pWCpJ2BTwEfBLYFRoDtgTdLuhv4VkQsrncASTsA1wGnR8Sm8QQZEfOB+QBDQ0Mxnn2YWY8YneBu9ErBE9w11VjNR9cCC4H3RMSz5RskHQKcIGnPiLis1g4kTSZLCN+NiO9XKbIOmF62vHtaZ2ZWnSe4K1TNpBARR42xbTmwfKwdSxJwGbA6Ii6sUWwRcKqkq4HDgI3uTzCzujzBXWFydTRL2h8YLC9f45d/uXcDJwAPSLo3rTsbmJE+fylwEzAbWAM8D5zUQOxmZtZkdZOCpAXA/sBK4JW0OoAxk0JE3AmoTpkg67cwM7MOkOdK4fA01sDMzHpcnttLl1YZdGZmZj0oz5XCQrLE8CTwElmTUETE/oVGZmZmLZcnKVxG6jDm1T4FMzPrQXmSwkhELCo8EjMza7s8SeEnkv4JuIGs+QjIdUuqmVk+S5d6MFqHyJMUXk+WDP64bF3dW1LNzHLxBHcdpW5SiAgPKDOz4lSb4M5JoW1q3pIq6fOSfmuM7UdI+kAxYZlZ3xid4G7SJE9w1wHGulJ4ALhB0ovACl6dJXVv4EDgR8BXCo/QzHqbJ7jrKMpmmhijgLQ32TxGU4EXgNXAkoh4ofjwtjY0NBTDw8PtOLSZWdeStDwihuqVy9On8FPgp02JyszMOlqRz2g2M7Mu46RgZmYlTgpmZlaS53kKA8BfsPVDdv68uLDMzKwd8oxo/hfgDrJbUF8uNhwz6wmetqJr5UkKb4iIzxUeiZn1Bk9b0dXy9CncKGl24ZGYWW+oNm2FdY08SeE0ssTwoqTN6bWp6MDMrEt52oqulmfw2o6tCMTMeoSnrehqefoUkHQM8N60eFtE3FhcSGbW9WbOdDLoUnWbjyR9lawJaVV6nSbp/KIDMzOz1svTpzAbOCoiFkTEAuBo4P31PiRpgaSnJD1YY/ssSRsl3Zte5zQWupmZNVuu5iNgF+CZ9H7nnJ+5HLgEWDhGmTsiws9kMDPrEHmSwvlkz2leDIisb+HMeh+KiCWSBicUnZmZtVSeu4+uknQb8K606nMR8WSTjj9T0n3Az4HPRMTKaoUkzQPmAcyYMaNJhzYzs0pjPY5zn/TvwWQP2FmbXruldRO1AtgjIg4ALgaur1UwIuZHxFBEDA0MDDTh0GZmVs1YVwqfJvt1fkGVbQEcMZEDR8Smsvc3SfqWpCkR8fRE9mtmTeZ5jPpKzaQQEfPS2/dFxIvl2yRtP9EDS3oL8IuICEmHkl21bJjofs2siTyPUd/Jc0vqXTnXvYakq4ClwDskrZV0sqRTJJ2SihwHPJj6FC4C5ka9B0abWWt5HqO+U/NKIf2Snwa8XtJBZHceAewEvKHejiPig3W2X0J2y6qZdarReYxGrxQ8j1HPG6tP4U+AE4HdgQvL1m8Gzi4wJjPrFJ7HqO+M1adwBXCFpGMj4roWxmRmncTzGPWVsZqPPhIR/wgMSvp05faIuLDKx8zMrIuN1Xz0xvTvDq0IxMzM2m+s5qPvpH+/1LpwzMysnfJMnf23knaSNFnSLZJGJH2kFcGZmVlr5Rmn8Mdp9PEHgMeAtwGfLTIoMyvQ0qVw/vnZv2YV8sySOlrm/cA/R8RGSWOVN7NO5RHKVkeeK4UbJT0EHALcImkAeLHOZ8ysE3mEstVRNylExJnA7wFDEfEb4FfAnKIDM7MCjI5QnjTJI5StqrrNR5ImAx8B3puajW4HLi04LjMrgkcoWx15+hS+DUwGvpWWT0jrPl5UUGZWII9QtjHkSQrvSg/CGXVrmtnUzMx6TJ6O5pcl7TW6IGlP4OXiQjIzs3bJc6XwWWCxpEfJps/eAzip0KjMzKwt6iaFiLhF0t7AO9KqhyPipWLDMjOzdshz99H2wCeB3yd7NvMdki6tfESnmbWJn6FsTZSn+Wgh2YN1Lk7LHwKuBI4vKigzy8kjlK3J8iSFd0bEvmXLiyWtKiogM2tAtRHKTgo2AXnuPloh6fDRBUmHAcPFhWRmuXmEsjVZniuFQ4C7JP0sLc8AHpb0ABARsX9h0ZnZ2DxC2ZosT1I4uvAozGz8PELZmijPLamPtyIQMzNrvzx9CuMiaYGkpyQ9WGO7JF0kaY2k+yUdXFQsZmaWT2FJAbicsZue3gfsnV7zyCbZMzOzNiosKUTEEuCZMYrMARZG5m5gF0lTi4rHrGv4cZnWRnk6mosyDXiibHltWre+sqCkeWRXE8yYMaMlwZm1hQejWZsV2XzUNBExPyKGImJoYGCg3eGYFcePy7Q2a2dSWAdML1vePa0z618ejGZt1s7mo0XAqZKuBg4DNkbEVk1HZn3Fg9GszQpLCpKuAmYBUyStBb5A9lhPIuJS4CZgNrAGeB4/o8Es48Fo1kaFJYWI+GCd7QF8qqjjm5lZ47qio9nMzFrDScGsaB53YF2knR3NZr3P4w6sy/hKwaxIHndgXcZJwaxIHndgXcbNR2ZF8rgD6zJOCmZF87gD6yJuPjIzsxInBTMzK3FSMGuUxx1YD3OfglkjPO7AepyvFMwa4XEH1uOcFMwa4XEH1uPcfGTWCI87sB7npGDWKI87sB7m5iMzMytxUjAzsxInBetvHnNg9hruU7D+5TEHZlvxlYL1L485MNuKk4L1L485MNuKm4+sf3nMgdlWCr1SkHS0pIclrZF0ZpXtJ0oakXRven28yHjMtjJzJpx1lhOCWVLYlYKkScA3gaOAtcA9khZFxKqKotdExKlFxWF9ZulS//I3m4Aim48OBdZExKMAkq4G5gCVScGsOXw3kdmEFdl8NA14omx5bVpX6VhJ90u6VtL0ajuSNE/SsKThkZGRImK1XuC7icwmrN13H90ADEbE/sDNwBXVCkXE/IgYioihgYGBlgZoXcR3E5lNWJHNR+uA8l/+u6d1JRGxoWzxH4C/LTAe63W+m8hswopMCvcAe0t6K1kymAt8qLyApKkRsT4tHgOsLjAe6weewdRsQgpLChGxRdKpwA+BScCCiFgp6VxgOCIWAX8l6RhgC/AMcGJR8ViX8t1EZi2liGh3DA0ZGhqK4eHhdodhreC7icyaRtLyiBiqV67dHc1mtfluIrOWc1KwzuW7icxaznMfWefy3URmLeekYJ3NdxOZtZSbj6y1/KQzs47mKwVrHd9NZNbxfKVgreO7icw6npOCtY7vJjLreG4+solpZMSx7yYy63hOCjZ+4+kj8N1EZh3NzUc2fu4jMOs5Tgo2fu4jMOs5bj6y8XMfgVnPcVKwV41nmmr3EZj1FCcFy3hgmZnhPgUb5U5jM8NJobc1Ms+QO43NDDcf9a5Gm4PcaWxmOCn0rmrNQR5YZmZ1uPmom7g5yMwK5iuFbuHmIDNrASeFdmpkXICbg8ysBZwU2qXRX/6jzUGj5d0cZGYFKLRPQdLRkh6WtEbSmVW2byfpmrR9maTBIuMpXCNt/o2OCxhtDjrvPA8sM7PCFHalIGkS8E3gKGAtcI+kRRGxqqzYycAvI+JtkuYCXwP+rJCAxjOFQyOfacUvfzcHmVnBimw+OhRYExGPAki6GpgDlCeFOcAX0/trgUskKSKiqZGMZwqHRj/TaJu/O4LNrAMV2Xw0DXiibHltWle1TERsATYCu1buSNI8ScOShkdGRhqPZDxTODT6mfHcAjpzJpx1lhOCmXWMruhojoj5wHyAoaGhxq8ixtNU0+hn/MvfzHpAkUlhHTC9bHn3tK5ambWStgF2BjY0PZLx/MEe72ecDMysixWZFO4B9pb0VrI//nOBD1WUWQR8DFgKHAfc2vT+hFHj+YPtP/Jm1mcKSwoRsUXSqcAPgUnAgohYKelcYDgiFgGXAVdKWgM8Q5Y4zMysTQrtU4iIm4CbKtadU/b+ReD4ImMwM7P8PCGemZmVOCmYmVmJk4KZmZU4KZiZWYmKugO0KJJGgMfH+fEpwNNNDKcbuM79wXXuDxOp8x4RMVCvUNclhYmQNBwRQ+2Oo5Vc5/7gOveHVtTZzUdmZlbipGBmZiX9lhTmtzuANnCd+4Pr3B8Kr3Nf9SmYmdnY+u1KwczMxuCkYGZmJT2ZFCQdLelhSWsknVll+3aSrknbl0kabH2UzZWjzidKGpF0b3p9vB1xNoukBZKekvRgje2SdFH6Pu6XdHCrY2y2HHWeJWlj2Tk+p1q5biJpuqTFklZJWinptCpleupc56xzcec6InrqRTZN9yPAnsC2wH3AvhVlPglcmt7PBa5pd9wtqPOJwCXtjrWJdX4vcDDwYI3ts4EfAAIOB5a1O+YW1HkWcGO742xynacCB6f3OwL/WeW/7Z461znrXNi57sUrhUOBNRHxaET8GrgamFNRZg5wRXp/LXCkJLUwxmbLU+eeEhFLyJ7BUcscYGFk7gZ2kTS1NdEVI0ede05ErI+IFen9ZmA1Wz/rvafOdc46F6YXk8I04Imy5bVs/YWWykTEFmAjsGtLoitGnjoDHJsur6+VNL3K9l6S9zvpNTMl3SfpB5L2a3cwzZSaeQ8CllVs6tlzPUadoaBz3YtJwaq7ARiMiP2Bm3n1Ssl6xwqy+W0OAC4Grm9zPE0jaQfgOuD0iNjU7nhaoU6dCzvXvZgU1gHlv4J3T+uqlpG0DbAzsKEl0RWjbp0jYkNEvJQW/wE4pEWxtUue/w56SkRsiojn0vubgMmSprQ5rAmTNJnsj+N3I+L7VYr03LmuV+ciz3UvJoV7gL0lvVXStmQdyYsqyiwCPpbeHwfcGqn3pkvVrXNFG+sxZO2UvWwR8NF0Z8rhwMaIWN/uoIok6S2jfWOSDiX7/7ubf+yQ6nMZsDoiLqxRrKfOdZ46F3muC31GcztExBZJpwI/JLsrZ0FErJR0LjAcEYvIvvArJa0h67ib276IJy5nnf9K0jHAFrI6n9i2gJtA0lVkd2BMkbQW+AIwGSAiLiV7NvhsYA3wPHBSeyJtnhx1Pg74hKQtwAvA3C7/sQPwbuAE4AFJ96Z1ZwMzoGfPdZ46F3auPc2FmZmV9GLzkZmZjZOTgpmZlTgpmJlZiZOCmZmVOCmYmVmJk4L1NUm7SPpk2fJukq4t6Fj/fazZLCX9rqTLizi2WV6+JdX6Wppb5saIeGcLjnUXcExEPD1GmR8Bfx4RPys6HrNqfKVg/e6rwF5pTvqvSxocfV5BegbF9ZJulvSYpFMlfVrSTyTdLem3Urm9JP2bpOWS7pC0T+VBJL0deGk0IUg6XtKDaUKzJWVFb6DLB1Nad3NSsH53JvBIRBwYEZ+tsv2dwP8A3gX8DfB8RBwELAU+msrMB/4yIg4BPgN8q8p+3k02idmoc4A/SROaHVO2fhh4zwTqYzYhPTfNhVmTLU5z2m+WtJHslzzAA8D+aSbL3wP+ueyRHNtV2c9UYKRs+T+AyyV9Dyif8OwpYLcmxm/WECcFs7G9VPb+lbLlV8j+/3kd8GxEHFhnPy+QzcYLQEScIukw4P3AckmHRMQGYPtU1qwt3Hxk/W4z2SMPxyXNc/9fko6H0vOCD6hSdDXwttEFSXtFxLKIOIfsCmJ06ue3A1WfwWzWCk4K1tfSr/P/SJ2+Xx/nbj4MnCzpPmAl1R+FugQ4qOyxr1+X9EDq1L6L7LnaAH8I/Os44zCbMN+SatYikv4euCEiflRj+3bA7cDvp8fEmrWcrxTMWucrwBvG2D4DONMJwdrJVwpmZlbiKwUzMytxUjAzsxInBTMzK3FSMDOzEicFMzMr+f81+KcIumYK6QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "m=1.0    # assume one kilogram\n",
    "F0=1.0   # and one newton, just for illustration\n",
    "\n",
    "v0 = 0.0  # start at rest\n",
    "x0 = 0.0  # at the origin\n",
    "xf = 3.0  # go to 3.0m\n",
    "dt = 0.1  # 0.1 sec intervals\n",
    "t = 0.0   # start at t=0.0s\n",
    "\n",
    "s=array([x0, v0])  # the \"state\" will be position and velocity\n",
    "\n",
    "def derivs_F(s, t):\n",
    "\n",
    "    x=s[0]     # extract position and velocity from the \"state\"\n",
    "    v=s[1]\n",
    "    \n",
    "    dxdt=v     # use the recipe here to get dvdt \n",
    "    dvdt=F0/m  # and dxdt\n",
    "    \n",
    "    return array([dxdt, dvdt])\n",
    "\n",
    "def HeunStep(s, t, derivs, dt):\n",
    "    f1=derivs(s,t)\n",
    "    f2=derivs(s+f1*dt,t+dt)\n",
    "    return s + 0.5*(f1+f2)*dt\n",
    "\n",
    "xlist = [x0]\n",
    "tlist = [t]\n",
    "\n",
    "while s[0] < xf:\n",
    "    s = HeunStep(s, t, derivs_F, dt)\n",
    "    t += dt\n",
    "    xlist.append(s[0])\n",
    "    tlist.append(t)\n",
    "\n",
    "title('Motion with constant force')\n",
    "xlabel('time (s)')\n",
    "ylabel('position (m)')\n",
    "plot(tlist, xlist, 'r.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.8.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}