{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Using a GAP model\n",
    "\n",
    "First fit a new potential using the `teach_sparse` program"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline\n",
    "import numpy as np\n",
    "import quippy\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "import os"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pot = quippy.Potential(\"\",param_filename=\"gp_new.xml\")\n",
    "potSW = quippy.Potential(\"IP SW\")\n",
    "potSW2 = quippy.Potential(\"IP SW\",param_filename=\"ip.parms.SW2.xml\")\n",
    "potSW3 = quippy.Potential(\"IP SW\",param_filename=\"ip.parms.SW3.xml\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7f8fa49a0110>]"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgoAAAFkCAYAAABB1xPiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3XmcXHWd7//Xt5d0QpZO0kknkJCGJBAhCKEbCIIgqAR0\nrjA4uDSo944goIgSQEZgYNS5gDBAcEGc++OHikALjDeyyJiw6QBDBNPshOwL2fdA9l7O/aO6QxK6\nqro7VXWqul5PHv14dJ/61jmffDlJvft8v+d8QxRFSJIkdaQk7gIkSVL+MihIkqSkDAqSJCkpg4Ik\nSUrKoCBJkpIyKEiSpKQMCpIkKSmDgiRJSsqgIEmSkjIoSJKkpLIaFEIIJSGEfw0hLAghbA0hzAsh\n/HM2jylJkjKnLMv7/z5wEfA14G3gGODXIYSNURT9PMvHliRJ+yjbQeFjwCNRFP2p7eclIYRzgeOy\nfFxJkpQB2Z6j8N/Ap0IIhwCEEI4CTgSeyPJxJUlSBmT7isKPgQHAOyGEFhLB5Nooin7XUeMQQhVw\nOrAI2J7l2iRJ6kl6AwcB06IoWpepnWY7KHwJOBf4Mok5ChOAn4QQlkdR9NsO2p8O3J/lmiRJ6snO\nAx7I1M6yHRRuAW6Koujhtp/fCiEcBFwNdBQUFgHweWBIYsOoGaOY+qup2a2yB5g8eTJTpkyJu4yC\nYp91j/3WdfZZ99hvXTNr1iy+8pWvQPtnaYZkOyjsB7Tsta2V5HMjEsMNQ4ADEht69+tNbW1tdqrr\nQSorK+2nLrLPusd+6zr7rHvst27L6NB9toPCY8A/hxCWAm8BtcBk4O4sH1eSJGVAtoPCt4F/Be4E\nqoHlwF1t2yRJUp7LalCIomgLcHnblyRJKjCu9dBD1NfXx11CwbHPusd+6zr7rHvst/wQoiiKu4Zd\nQgi1wEwuZNdkxkOnH8rsF2bHWZYkSXmvsbGRuro6gLooihoztd9sz1HoloF/GUh1VTUANdU1MVcj\nSVLxysugcNKFJ/Ho5Y/GXYYkSUUvL+cozF8/P+4SJEkSeRoUFmxYQGvUGncZkiQVvbwMCtubt7Nk\n05K4y5AkqejlZVAAeGv1W3GXIElS0cvLoLBf+X68ufrNuMuQJKno5WVQGD1oNG+t8YqCJElxMyhI\nkqSk8jIojBk8hllrZnnngyRJMcvPoDBoDNuat7Fww8K4S5EkqajlbVAAHH6QJClmeRkUhvYdSmVF\npbdISpIUs7wMCiEExleP94qCJEkxy8ugADB+qEFBkqS45XVQeGftO7S0tsRdiiRJRSt/g0L1eLY3\nb2fBhgVxlyJJUtHK36AwdDzgnQ+SJMUpb4PC8H7DGdR7kHc+SJIUo7wNCt75IElS/PI2KEBi+MFV\nJCVJik/eB4XZ62bT3NocdymSJBWl/A4K1ePZ2bKTeevnxV2KJElFKb+DQvudD05olCQpFnkdFKr7\nVjNkvyFOaJQkKSZ5HRRCCD7KWZKkGOV1UIC2NR8cepAkKRb5HxSqxzNn3RyaWpriLkWSpKKT/0Fh\n6HiaWpuYu35u3KVIklR08j8oVHvngyRJccn7oDBkvyFU9612QqMkSTHI+6AAeOeDJEkxKZyg4NCD\nJEk5VxhBoXo8c9fPZWfLzrhLkSSpqBRGUBg6nubWZuasmxN3KZIkFZXCCAptdz645LQkSblVEEFh\ncJ/BDO833HkKkiTlWEEEBfDOB0mS4lAwQeGI6iMMCpIk5VjWg0II4YAQwm9DCGtDCFtDCK+FEGq7\nup/xQ8czb/08tjdvz0aZkiSpA1kNCiGEgcALwA7gdOAw4ApgQ1f3Nb56PK1RK7PXzs5skZIkKamy\nLO//+8CSKIou2G3b4u7s6PChhwPw1pq3OGr4URkoTZIkpZPtoYfPAX8LITwUQlgVQmgMIVyQ9l0d\nGNh7ICP6j/DOB0mScijbQWE08E1gNjAJuAv4aQjhq93Z2fhq73yQJCmXsj30UAK8FEXRdW0/vxZC\nOAK4GPhtsjdddtllDBw4cI9t9fX1jB86nsfmPJa1YiVJKgQNDQ00NDTssW3Tpk1ZOVa2g8IKYNZe\n22YBn0/1posuuojzzjvvQ9u3Nm7ljhl3sK1pG33K+2SuSkmSCkh9fT319fV7bGtsbKSuri7jx8r2\n0MMLwLi9to0jzYTGhx9+ssPt46vHExHxztp3MlOdJElKKdtBYQpwfAjh6hDCmBDCucAFwM9Tvem5\n5xYzduwZVFeP4/rrb9q1ffc7HyRJUvZlNShEUfQ34GygHngDuBb4bhRFv0v1vvXrb2Pz5nKuu+4q\nrrvuyl3bB1QM4MABB3rngyRJOZLtOQpEUfQE8ERX31dZGbj00vM/tH189XjeXOMqkpIk5UJervVQ\nUjKV1tbmDl8bP3S8VxQkScqRvAwKAwYEamsP7vC18UPHs3DjQrbs3JLjqiRJKj55GRQ2bvx7brzx\nzg5fO6L6CABmrd37rktJkpRpeRkUSkth2rSOXzts6GEADj9IkpQDWZ/M2B0TJsCf/gTf+tYH2yad\nPYnFqxOPXyjbUMbkqZO5se+Nu16vqa5h+tTpuS5VkqQeLS+DwgknwD33wI4dUFGR2LZ49WLmTJqz\nq82Gtv92MSNIkpRxeTn0cMIJsGULvPBC3JVIklTc8jIoHHIIDB+eGH6QJEnxycugEAKccYZBQZKk\nuOVlUIBEUHjjDVi2LO5KJEkqXnkbFD79aSgpSX6bpCRJyr68DQpVVXDccQ4/SJIUp7y8PbLdGWfA\nHXdAc3PiOQm73wK5ecdmlr+/nNGDRlNWWpZ4XZIkZVTeB4Uf/ABeeokPPUxp5eaV7H/b/tz8hZs5\n5/Bz4ilQkqQeLm+HHgCOOQYGD+54+GF4v+GMqhzFX5f+NfeFSZJUJPI6KJSWwqRJyecpTBwxkRnL\nZuS2KEmSikheBwVIDD/87W+wZs2HX5s4YiIzl8+kqaUp94VJklQE8j4oTJoEUQTTO1jL4fiRx7Ot\neRtvrn4z94VJklQE8j4o7L//B6tJ7q12/1rKSsr46zLnKUiSlA15HxQgMfwwbRq0tu65vU95H44c\ndqRBQZKkLCmIoPCZzyTmKLzyyodfmzhiIjOWOqFRkqRsKIig8LGPQf/+HQ8/TBwxkXfWvsPG7Rtz\nX5gkST1cQQSF8vLE2g8dBoWREwF4ednLOa5KkqSeryCCAiTmKbz4Imzc68LBoVWHMrD3QOcpSJKU\nBQUTFE4/HVpa4Omn99xeEko4bsRxBgVJkrKgYIJCTQ0cdljyeQozls4giqLcFyZJUg9WMEEBEsMP\nf/pT4gFMu5s4YiJrt65l4caF8RQmSVIPVXBBYelSePvtPbcfN+I4ABeIkiQpwwoqKJx8MvTp8+Hh\nh6F9hzJm0BjnKUiSlGEFFRR694ZTTkl+m6RBQZKkzCqooACJ4Ye//CXikkuu2WPy4sQRE2lc0ciO\n5h0xVidJUs9SkEGhqWkm//7vP6OxsXHX9okjJrKzZSevrXotxuokSepZCi4oHHII9Ov3MC0tt3HX\nXQ/v2j5h+AR6lfZyQqMkSRlUMEHh+utvorp6HIce+hlaWuYA3+Dxx2czduwZVFeP44Yf3c6E4ROc\npyBJUgaVxV1AZ1133ZUMHVrNDTf8gW3bpgKwatVU4HNcd91VXHzx19j01EqemPdEvIVKktSDFMwV\nhfLyci699HwqK8Me2ysrA5deej7l5eVMHDmReevnsXbr2piqlCSpZymYoNCutbWZPn3upqLiNEK4\nm9bW5l2vTRyRWEnypWUvxVWeJEk9SsEFhdrag7n55sBttz1CFAXGjz9412ujB41myH5DnNAoSVKG\nFMwchXYPPngnACtXwqWXns/nP3/+rtdCCK4kKUlSBhXcFYV2w4fDxInw6KN7bp84YiIvLXuJ1qg1\nnsIkSepBchYUQgjfDyG0hhBuz9Q+zzwz8Tjn7ds/2Hb8yOPZsH0Dc9fNzdRhJEkqWjkJCiGEY4EL\ngYw+NvHMM2HLFnj22Q+27VpJ0uEHSZL2WdaDQgihH3AfcAGwMZP7PvxwGD16z+GHgb0HMq5qnBMa\nJUnKgFxcUbgTeCyKomcyveMQElcVHn0UdlsfypUkJUnKkKwGhRDCl4EJwNXZOsZZZ8Hy5bDb+lBM\nHDGR11a9xrambdk6rCRJRSFrQSGEMBK4AzgviqKmbB3nxBNh0KA9hx+OH3k8za3NvLLylWwdVpKk\nopDN5yjUAUOBxhBC+3OXS4GTQwjfBiqiaPcBgw9MnjyZysrKPbbV19dTX1//obbl5fDZz8Ijj8AP\nf5jY9tHqj9K7rDczls7ghANPyNgfSJKkfNDQ0EBDQ8Me2zZt2pSVY4Ukn9X7vuMQ+gI1e23+NTAL\n+HEURbM6eE8tMHPmzJnU1tZ2+lgPPQRf+hIsWgQ1bUf8+D0fZ8SAETx4zoPd+wNIklRAGhsbqaur\nA6iLoqgxXfvOytrQQxRFW6Ioenv3L2ALsK6jkLAvTj89cWXhscc+2DZxxETvfJAkaR/l+smMWbl8\nUVkJp5yy5zyFiSMnsnjTYlZtXpWNQ0qSVBRyGhSiKPpkFEWXZ2PfZ50Ff/4ztA/RtK8k6W2SkiR1\nX8Gu9bC3z30Omppg2rTEz6MqRzG833BmLJ0Rb2GSJBWwgls9MplRo2DChMTdD3c3TGLx6sW8/977\n3BHdwe8rf79H25rqGqZPnR5TpZIkFY4eExQg8ZTGn/4Uhh62mLmnz9m1fQ5z9mxoRpAkqVN6zNAD\nJILCxo2wzQcySpKUET0qKNTWwogRsHlz3JVIktQz9Kig0L5IlEFBkqTM6FFBARJBoTlrK0tIklRc\nelxQOPVUCD3uTyVJUjx63EdqRQX03S/uKiRJ6hl61O2R7WqG1/DWz2D0GCgrg2XvLaM1auXAygMT\nr1fvvVaVJEnqSI8MCn/543Sqq+Hqb8AFF8BdL9/Fd/70HV6+6mUGVAyIuzxJkgpGjxt6AKiqgo9/\n/INFoiaNmURzazPPLnw23sIkSSowPTIoQOLuhyefhC1bYMzgMYwZNIZp86fFXZYkSQWlRweF7dvh\nyScjLrvsGk4bcxrT5/vsZkmSuqLHBoVDDoGPfAR+9auZ3HnnzxjLWOZvmM/89fPjLk2SpILRY4MC\nwFlnwbRpD9PcfBuv/WEpZSVlXlWQJKkLemRQuP76m6iuHsd9932GHTvmAN9g+mOLKFvRn8t+/j2u\nv/6muEuUJKkg9MigcN11V3LddVfR3FwGTAUCq1ZNpXRRNWF0C9+/5rtxlyhJUkHokUGhvLycSy89\nn8rKsMf2qg3V7GA7M1fNjKkySZIKS48MCu1aW5upqLgbOI2KirspX9eXwX0GO09BkqRO6tFBobb2\nYG65JTBixCMcd1yg7ujRnDb6NJ+nIElSJ/XIRzi3e/DBOwF49124997zeeaZ87n3jXu44NELWLd1\nHVX7VcVcoSRJ+a1HX1Fo9+Uvw+rV8Oc/Jx7nHBHx1IKn4i5LkqS8VxRBobYWxo6F3/0ORg4YyeFD\nD3eegiRJnVAUQSGExFWF3/8edu6E08eczrT504iiKO7SJEnKa0URFCARFDZuhOnTE8MPy95fxqy1\ns+IuS5KkvFY0QWH8eDjiiMTww8k1J1NRWsG0ed79IElSKkUTFCBxVeGRR4Cm/Tip5iSmL3CegiRJ\nqRRVUPjSl2DzZnjiicQ8hb8s+gvbm7fHXZYkSXmrqILC2LFwzDGJ4YfTx5zOtuZtPL/k+bjLkiQp\nbxVVUIDE8MMf/wijeh/B/v32d56CJEkpFF1Q+OIXYft2eOyxwKQxk5ynIElSCkUXFA48ED7+8cTw\nw6Qxk3h91euseH9F3GVJkpSXii4oQGL4Ydo0OGbwaQA8ueDJmCuSJCk/FWVQOOccaG2F//rPodTu\nX+tqkpIkJVGUQWHYMPjkJz+4++HJ+U/SGrXGXZYkSXmnKIMCJIYfnn0Wjhk0iTVb1/DqylfjLkmS\npLxTtEHh7LOhtBSWvHACfcv7upqkJEkdKIu7gLgMHgyVIyZx9S2LKR0IP3rgR/yq8ld7tKmprmH6\nVAOEJKl4FW1QACjrt5jt58zZ9fMc5uzZwIwgSSpyRTv0ANCvb9wVSJKU37IaFEIIV4cQXgohvBdC\nWBVCmBpCODSbx+yKktK4K5AkKb9l+4rCScDPgInAp4FyYHoIoU+WjytJkjIgq3MUoij67O4/hxD+\nF7AaqANctlGSpDyX6zkKA4EIWJ/j40qSpG7IWVAIIQTgDuD5KIreztVxJUlS9+Xy9shfAIcDJ6Zr\nOHnyZCorK/fYVl9fT319fUYLanq/ldJflFNSUk5TU4DSITBwMaVbK2BnK03DfKyzJCn/NDQ00NDQ\nsMe2TZs2ZeVYIYqirOx4j4OE8HPgc8BJURQtSdGuFpg5c+ZMamtrs15XU1MTv/zlvdxwwx9Yteox\nIIJvHUHvjRu55fgfcfHFX6O8vDzrdUiStK8aGxupq6sDqIuiqDFT+8360ENbSDgLODVVSIhDeXk5\nl156PpWVoW1LgDfOY8foVXz9oi8bEiRJRS/bz1H4BXAecC6wJYQwrO2rdzaP21Wtrc306XM3ZWWn\nEd4KRGUtPDr70bjLkiQpdtm+onAxMAD4M7B8t68vZvm4XVJbezA33xy47bZHYEM1gzYP44E3H4i7\nLEmSYpft5ygUxCOiH3zwTgDeew+uueZ8Tui/jT/Nm8y6reuo2q8q5uokSYpPQXyQ58qAAVBfD6/c\n9wWiKOLhtx+OuyRJkmJlUNjLhRfC8rnDOKr/p3ngDYcfJEnFzaCwl2OOgaOPhuiNc3luyXMs2ZRX\nN2pIkpRTBoW9hJC4qvDa786morQ3v3vzd3GXJElSbAwKHTj3XOhT0p9DWs50+EGSVNQMCh1on9S4\n8slzeW3Va7y1+q24S5IkKRYGhSQuugjWzjiDvqUDvaogSSpaBoUk6urg6CMrqFp1Dg+8+QC5WBND\nkqR8Y1BIIoTEVYV3nziXRRsXMWPpjLhLkiQp5wwKKdTXQ5/VJ9OfAxx+kCQVJYNCCgMGwHn1pfBG\nPQ++9SDNrc1xlyRJUk4ZFNK48EJ4/7/PZc3WNTy94Om4y5EkKacMCmnU1cGE4UfTb/s4V5SUJBUd\ng0IaIcDFFwW2zDiX37/9f9natDXukiRJyhmDQifU10PvufVsadrM43Mej7scSZJyxqDQCQMGwFc+\newi91hzL/a87/CBJKh5lcRdQCCadPYnZSxazc9UGHn3yZcb+dCylJaW7Xq+prmH61OkxVihJUnYY\nFDph8erFLDlzzq6f5zN/zwZmBElSD+XQgyRJSsqgIEmSkjIodMLatev36XVJkgqVQaETqqoG7dPr\nkiQVKoNCJ4QQ9ul1SZIKlUFBkiQl5e2RnVBTXQPTYfnyVfTp05uBlQNY+O56mvuuY3i/4dSMqom7\nREmSssKg0AkdPUzpqafgtF/9HX0mvMu0K6fFUJUkSdnn0EM3fepT8JH1V7Jw6xs8ueDJuMuRJCkr\nDArdFAL8+OJTYMXR/PMTt8VdjiRJWWFQ2Aef+1zgwHev5OX103l91etxlyNJUsYZFPZBSQnc+JUv\nwKaRXPP47XGXI0lSxhkU9lH9F8sZMu8ynljyAMvfXx53OZIkZZRBYR+VlsK/nHkBUVNvrv/jz+Iu\nR5KkjDIoZMCFX6uk/5wLufftX7J55+a4y5EkKWMMChnQqxd87+Tv0BTe5+bp98RdjiRJGWNQyJAr\nvzGK3vO+xB1/nUJza3Pc5UiSlBEGhQzp0wcuPPIKNpct4v97fmrc5UiSlBEGhQz639+qpezdU/nh\nk7cSRVHc5UiStM8MChnUvz988cArWFX2Eo+99kLc5UiStM9cFCrDlr18B7zTi3MeOIODh4z40Os1\n1TUdLjIlSVI+Mihk2IoNS+CrO2liJ3OY8+EGZgRJUgFx6EGSJCWV9aAQQrgkhLAwhLAthDAjhHBs\nto8pSZIyI6tDDyGELwG3ARcCLwGTgWkhhEOjKFqbzWPHZe3a9fv0uqTCN+nsSSxevTjxQwSr16yl\neugQCIlN7XOVOtsuG/ss1mOr67I9R2Ey8O9RFN0LEEK4GPg74OvALVk+diyqqgaxnuQZqKpqUA6r\nkZRJnf0wWrx6MXMm7TlHaSO7/ZLQ9nnV2XZdaeux0+9TXZO1oBBCKAfqgBvbt0VRFIUQngI+lq3j\nxi2EsE+vS8qtrvwm6oeRilE2rygMAUqBVXttXwWMy+JxJanTMvnhv31nM4+8+CZLV+79z96eFqxY\nyvFX/wvzV7ybst28FUsY/71v7/o+lfkrl1B3/SXMX5m+3bE/6Fy7CddeQhSRtu28FUsYf9Ulnarx\nyKs7d+y66y7Z9X26tkddk36fDvt2X17eHjl58mQqKyv32FZfX099fX1MFXVeTXUNTIfly1fRp09v\nKiv7M2/NYqCFsUPHUDO8Ju4SpaKQqTHr5pYWHpv5Mo0LFrFo5dKUbZe8t4C/n/7RtP+yNpftYGbr\nPbSUNaVsF5W2sqj1eYgCUWlryrYtoZU3N75IS0jf7rV1nWhX0so7m18khMT3aetsfjF9jSWtzN32\nYtr9tYRW3tz04q7v0+1zztb0++xpw74NDQ00NDTssW3Tpk1ZOVY2g8JaoAUYttf2YcDKVG+cMmUK\ntbW12aorqzr6h+c/nnudLzx1NB+puJjHrrk8hqqk4tPZKwXpftNcsGE+Zz5+XOKH0tRDh1XlI7jx\n2P/gh4/Us5xFSdsdOnQMs2+ezbjnx3X8vJU2h1QfxOzbXgVg3F9Ttz102EHM/mkj42Z2ot3PGhl3\nYpp21Qcx+yeNiWO/3Ik6b29MX2P1Qcy+ozHt/tr/LED6P08n99nThn07+uW5sbGRurq6jB8ra7dH\nRlHUBMwEPtW+LST+T30K+O9sHTcfnXPSkXx058U8vvmHvDo39SVJSbmx9r33OOmGy9gYbU7ZbkDJ\n/vx0/Exe/Pt1HFJ1SMq2VZV9ufCzx9OvT69MlirFKttDD7cDvw4hzOSD2yP3A36d5ePmnT9850eM\n/WkD//CLa5k/5e64y5EKVmeHFNJdKVi/cyUvrHmc0FoObE/abvjg/lx6TuIKZ/ARdSpCWQ0KURQ9\nFEIYAvyIxJDDq8DpURStyeZx89Ho/as474B/5b51l/Lvj32Tiz6X+ctDUjHozJDC60sW01ya+lN9\nZN/RLLx1HuNfHscc3s9oje1zlVK+3oV22dhnsR5bXRfyaTnkEEItMHPmzJkFO0chlZ3NzQz6fi3s\n6M+G256nV6+eNWYm5cK4E8d9KCjsLtzfi+i8nfAAcG7y/Rw6/VBmvzA77f7a28FeVzM64EN9FKfd\n5ijURVHUmKn95uVdDz1Vr7Iybp/0Ey5+8ZOcf8cD/Paq8+IuSSo4aW9zaw58Y+DDTB9wNYuZl3Z/\nXflN1BCgYmRQyLGLJp3KLc/+A/dvvoofvHsWYw7sF3dJUuw6O+9g+84mQq/UEwXHDh/F//nuOYx7\n6NpOHdsPfyk1g0IM/uPCW6m95zDOvv0mXp9yQ9zlSLFLN+/gvanbOfqa7/J6awOtO1NPcWq/Dc4x\naykzDAoxOPrgg/i7gd/jj9zCw0+dzxc+PTrukqS8tnLbElY3P8RRJV9lRZ/fs5Lk8wTaeaVAygyD\nQky2/tfzMKuZLz9wJNcOG8EaVzpTEUs376BP8wA2/O+lVPQqZdxzj6d+YpukjDIoxGTZ2mVwXgut\nbGFu29PEXFxGPU1n5x706dsPUqy6OnLYMCp6le56j0MKUu4YFCRlTbq5B2sf3sLQ75zJ2i2LUu5n\n98fveqVNyi2fMxaTdJdaXelMxWB90zK2li9hQOneS8JIyhcGhZikW8msp610puKULvD2balk862v\nMHxQZcp2kuLj0ENM0q1k1tNWOlNxqho8iPUp5h4cUF1NCMF5B1IeMyhI6rJ0kxQPqDqQAz/xNeau\nTX0bY3sgdt6BlL8MCvkqf5bgkD4k3STFOffPh7qnKQ19aWFnrsuTlEEGhZgku9S6dMMqtrZsoqm/\ncxRUuHrRl8fPepFvTz+bOSRfcElS/jMoxCTZpdYdzTsY+S8nsnDLGp57eSMnHTswx5VJ6aWbpNgv\n9OK0CYc790DqAQwKeaairILnLn2Y8T8/ms/88h9Zesj/ZeBAJzYqf0QRlFb0Sdmm/a4d5x5Ihc+g\nkIc+Mvxg7vzUb/jmf/09n/j+Hbx612S8CULZlm6C4qjqUZx8zj9x619v5L0d76bcl3ftSD2HQSFP\nXXzqWfzn21fyaMtVXDHleG6//GNxl6QeLt0Exbn3LeapCU/Rr7yOwb0OYD3Lc12ipBj4wKU89h8X\n38j+0XFMWfZFpj+f/F50KRcCJdw2YTrv3foyQwb0i7scSTniFYU8Vl5azuhXSlkxZzmn/6GGMUNG\nsG79OleZVFakm6A4sKQvl591GuDCTFIxMSjkuTUbVsF5rcBW5jMXcJVJZUflwMqUT1Hc/bHihlOp\neBgUpB4u3STFoZUHUHbsqSzcsDDlfpygKBUng0Kec5VJ7av0T1GcC7Uz6MUAdrIx1+VJynNOZsxz\nrjKpbOsdBvLGBYs4aGh13KVIykNeUchzrjKpfZXuqtN+lHLEQcOcoCipQwaFArfhPRfcUWq9+uyX\n8nWfoigpFYNCgVuz412m3DeLyV85LO5SlGMpJylGQEl/1h1bxbptS1Lux6tSklIxKOS5VJeDm1ta\nWFyxkstfOwV4lslfOTyXpSlm6SYp8gD05igGlu3PRlbkuDpJPYVBIc+luxy8YtMaDr/501z++ilE\n9z3D5V85IkeVKd9V9RrB6tte4bCTPmJQkNRtBoUCt3/lUGZ9/2kOu+k0rvhxLT+++QAGDajo8H55\nn+LYs6SbpBjt3EFJSXCSoqR9YlDoAYYPGMI7Vz/FAb8fwZrzFrNmt9d8imNhSfdwpPawN3PuUrZE\nLSn35SRFSZlgUOghhg2o4uDBI5nP/LhL0T5IN+/g/T/sYMz3vsqCPr8DWlPuy0mKkjLBBy71IBs2\nbEr5uk9xLHwrti5mcfQcZ/b5N0YPGh13OZKKgFcUepCqqkGdXtRH+SldmKto7s+GG+bRp6KMcY/f\nlaOqJBVXKlPzAAALvElEQVQzg0IP4lMcC1+6sDdq2HD6VCT+2jpJUVIuGBSKyIp1m2lthRIHnHIu\n3STFUUNHcUr9FSxYszzlfnYPe05SlJQLBoUi8n60nIMu+SZ/ufY2Dh6ZeKxvFEVMnnwtU6bc4BWH\nLEq7guN9C3nqnaeAXjmuTJJSMyj0ICkvRUewbVAJ71b9htGTfsuwsiFU9q9g+/btLFnyLlOfu5/e\nvXvv2o+/reZWGb35Sd2z/GTa15nDnPRvkKQcMSj0IJ35cP/LrLc59aGjWfWlxazabfsSdlsPwIyQ\ncekmKQ4oqeBb/+NE/vD/O+9AUn4xKBSZTxx2OANDfzawLmkbb6PsvHRzD0YOGcWJ50xmU+v2lPvx\n4UiS8pVBoQgNqRqcMih4G2XnpZ17cP98nlmQfu6B80Mk5auszH8PIdSEEO4OISwIIWwNIcwNIfwg\nhFCejeOpa9J9KDU3RzmqpOfrRX8aTpnJoUMPirsUSeqWbF1R+AiJi6/fAOYDRwB3A/sBV2XpmMqQ\nhRsXcNK1NxC98jRrNi1LbCyiRaY6u94CwMrVa5LsJaFfKOPLn6jlX7xgIKlAZSUoRFE0DZi226ZF\nIYRbgYsxKOS9ijCA50t+BPOa4bw91xMohkWm0g0nbHmkidN/9G+8sG4qW0o2pNxX+zCOD0eSVKhy\nOUdhIOAsuQJQM6Sah89/muMe+ig72Bh3ORnTlSsFqSzbspBlO69nRNkZlJUNZxMrk7ZtH+bpiVde\nJBWHnASFEMJY4NvA5bk4nlLrzG+3Rx40kr6UsSPFfna/OyJTH8LZlO5KQXufpLvro3dzf5ZetYKq\nAX0ZN2NcyqAgSYWuS0EhhHAT8E8pmkTAYVEU7frXOIQwAvhP4MEoiu7pVpXKqM5+YKdbd2BT63a+\neWcD//T5z3b6QzgbMhVSVm14n3GXf4v1ralX4Txw2HCqBvTd17IlqSB09YrCrcCv0rRZ0P5NCOEA\n4Bng+SiKLursQSZPnkxlZeUe2+rr66mvr+9CqdpX6e6OaKWJX649l1/eVU7J6s7d0NLZD/WufPhn\n6krBppYVbC95hl5hP3aSPCzs3i/OPZAUh4aGBhoaGvbYtmlT6l9yuqtLQSGKonWQ4gb83bRdSXgG\neBn4eleOM2XKFGpra7vyFsXgkOoa7qt/kn977A883HpFyrbtH9Kd/VDP5BWK9e9t5cwbfsb7rTtT\ntqsZcDCLbn2HcS+OY06KoLC7uIdTJBWnjn55bmxspK6uLuPHytZzFA4A/gwsJnGXQ3UIYVgIYVg2\njqf4HHvoKB664jscMuzglO3Wt25kyHc/w+I1qW8n7LQI1m3cyV/fWsHyVatTNl27cymPbb+SZlI/\nHbGil4/5kKS9ZWsy42nA6Lavd9u2BRJzGEqzdExlWFcuq6cbpuhd0p/y0l7siFL/pj5n1Xz6XjaR\nrSsXpG63dg5DflKR+CHNGTWiXw0LfzCPI54Z36kFlxxOkKQPZOs5Cr8BfpONfSt3MnlZfdSQocy+\n/REOnXEoc5mbtF2vkn6M6v1R5pXOopn3k7brWzKErw7/N0YOHsLP97uElbsvarV3294VlJd2/lR3\nOEGSPuBaD8qpdFceDho6jFk/vptxzz3HnBRBYUTVYO666H8BcO8dV3TqBkWvFEhS1xkUlBGF8CHs\nlQJJ6jqDgjIizg/hQggpklSoDArKqc5+qHflw98rBZKUPQYF5VRnP9T98Jek/JCV5yhIkqSewaAg\nSZKSMihIkqSkDAqSJCkpg4IkSUrKoCBJkpIyKEiSpKQMCpIkKSmDgiRJSsqgIEmSkjIoSJKkpAwK\nkiQpKYOCJElKyqAgSZKSMihIkqSkDAqSJCkpg4IkSUrKoCBJkpIyKEiSpKQMCpIkKSmDgiRJSsqg\nIEmSkjIoSJKkpAwKkiQpKYOCJElKyqAgSZKSMihIkqSkDAqSJCkpg4IkSUrKoCBJkpIyKEiSpKQM\nCpIkKSmDgiRJSsqgIEmSkjIoSJKkpAwKPURDQ0PcJRQc+6x77Leus8+6x37LD1kPCiGEXiGEV0MI\nrSGEI7N9vGLlX6ius8+6x37rOvuse+y3/JCLKwq3AEuBKAfHkiRJGZTVoBBC+AxwGnAlELJ5LEmS\nlHll2dpxCGEY8H+AM4Ft2TqOJEnKnqwFBeBXwC+iKHolhFDTyff0Bpg1a1b2quqhNm3aRGNjY9xl\nFBT7rHvst66zz7rHfuua3T47e2dyvyGKOj91IIRwE/BPKZpEwGHAGcA5wClRFLWGEA4CFgAToih6\nPcX+zwXu73RBkiRpb+dFUfRApnbW1aBQBVSlabYQeAj4H3ttLwWagfujKPrHFPs/HVgEbO90YZIk\nqTdwEDAtiqJ1mdppl4JCp3cawkhgwG6bDgCmAf8AvBRF0fKMH1SSJGVcVuYoRFG0dPefQwhbSNz1\nsMCQIElS4cjlkxl9joIkSQUmK0MPkiSpZ3CtB0mSlJRBQZIkJZXToBBCOCmE8GgIYVnbIlFnpmn/\nibZ2u3+1hBCqc1VznEIIV4cQXgohvBdCWBVCmBpCOLQT7zslhDAzhLA9hDAnhPA/c1FvvuhOvxX7\nuQYQQrg4hPBaCGFT29d/hxDOSPOeYj/XutRnnmcfFkL4fls/3J6mXVGfa3vrTL9l6nzL9RWFvsCr\nwLfo/OTGCDgEGN72tX8URauzU17eOQn4GTAR+DRQDkwPIfRJ9oa2h1s9DjwNHAX8BLg7hHBatovN\nI13utzbFfK4BvEvigWq1QB3wDPBICOGwjhp7rgFd7LM2xX6e7RJCOBa4EHgtTbuD8FzbpbP91mbf\nz7coimL5AlqBM9O0+QTQAgyIq858+gKGtPXbx1O0uRl4fa9tDcATcdef5/3mudZxv6wD/jHJa55r\nXe8zz7MP+qIfMBv4JPAscHuKtp5r3eu3jJxvhTBHIQCvhhCWhxCmhxBOiLugGA0kkQ7Xp2hzPPDU\nXtumAR/LVlEFoDP9Bp5ru4QQSkIIXwb2A15M0sxzbTed7DPwPGt3J/BYFEXPdKKt59oHutJvkIHz\nLZuLQmXCCuAi4G9ABfAN4M8hhOOiKHo11spyLIQQgDuA56MoejtF0+HAqr22rQIGhBAqoijaka0a\n81EX+s1zDQghHEHiQ6438D5wdhRF7yRp7rlGl/vM8wxoC1QTgGM6+RbPNbrVbxk53/I6KERRNAeY\ns9umGSGEMcBkoNgmsvwCOBw4Me5CCkyn+s1zbZd3SIwBV5JY2O3eEMLJKT741IU+8zzb9Yj/O4BP\nR1HUFHc9haI7/Zap860Qhh729hIwNu4icimE8HPgsyRW41yRpvlKYNhe24YB7xVL6m7XxX7rSNGd\na1EUNUdRtCCKoleiKLqWxGSp7yZp7rlGl/usI8V2ntUBQ4HGEEJTCKGJxFj6d0MIO9uuAu7Nc617\n/daRLp9veX1FIYkJJC6nFIW2D7uzgE9EUbSkE295EfjMXtsmkXrMtMfpRr91pKjOtSRKSFyy7Ijn\nWsdS9VlHiu08ewr46F7bfg3MAn4ctc3C24vnWvf6rSNdPt9yGhRCCH1JJJn25DM6hHAUsD6KondD\nCDcBB0RR9D/b2n+XxLLVb5EY//sGcCpQFLfEhBB+AdQDZwJbQgjtiXpTFEXb29rcCIxo7zPgl8Al\nIYSbgXuAT5G4HPrZnBYfo+70W7Gfa7CrT/4TWAL0B84j8RvLpLbX9/j7iedal/vM8wyiKNoC7DFf\nKCQWDlwXRdGstp/9d20v3em3TJ1vub6icAyJ2zmitq/b2rb/Bvg6iQkrB+7WvldbmwOArcDrwKei\nKPqvXBUcs4tJ9NOf99r+j8C9bd/vz259FkXRohDC3wFTgO8AS4Hzoyjae8ZwT9blfsNzDaCaxN/F\n/YFNJPpg0m6zq/f4++m5BnSxz/A8S2bv34b9d61zUvYbGTrfXBRKkiQlVYiTGSVJUo4YFCRJUlIG\nBUmSlJRBQZIkJWVQkCRJSRkUJElSUgYFSZKUlEFBkiQlZVCQJElJGRQkSVJSBgVJkpTU/wNDF5V4\nj5HlsgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f8fbdbd0310>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lattice = np.eye(3)*20.0\n",
    "at = quippy.Atoms(n=2,lattice=lattice)\n",
    "at.set_atoms([14 for i in range(at.n)])\n",
    "at.pos[:,1] = [0.0, 0.0, 0.0]\n",
    "\n",
    "GAP2bEnergy = []\n",
    "SW2bEnergy = []\n",
    "rs = np.linspace(1.5,4.1,40)\n",
    "for r in rs:\n",
    "    at.pos[:,2] = [r,0.0,0.0]\n",
    "    pot.calc(at,energy=True,force=True)\n",
    "    GAP2bEnergy.append(at.energy)\n",
    "    potSW.calc(at,energy=True,force=True)\n",
    "    SW2bEnergy.append(at.energy)\n",
    "    \n",
    "plt.plot(rs,GAP2bEnergy,\"*-\")\n",
    "plt.plot(rs,SW2bEnergy,\"s-\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pot.calc(at,energy=True,args_str=\"energy_per_coordinate\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7f8fa4665b90>]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgsAAAFkCAYAAACuFXjcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3Xl8VOXd///XJwsmgAQkhN0AgmjtXTXRKri0bhHrigsa\ntVWk3lKVaqzVX2+/6Fd7K1+1SmsFobdWsbWp2hZFbBW1m7s2Ue+2KqBAWAIJIRC2JGS5fn/MIElM\nTmYmkzmzvJ8+8niQc64z5zPHA/Oec13nOuacQ0RERKQraX4XICIiIvFNYUFEREQ8KSyIiIiIJ4UF\nERER8aSwICIiIp4UFkRERMSTwoKIiIh4UlgQERERTwoLIiIi4klhQURERDyFHRbM7AQzW2JmG8ys\n1czOCWGbPmZ2t5mtMbMGM1tlZldGVLGIiIjEVEYE2/QDPgQeA/4Q4jbPAkOA6cDnwHB0VUNERCQh\nhB0WnHMvAS8BmJl1197MpgAnAOOcc9uCi9eGu18RERHxRyy+3Z8N/AO41czWm9lyM7vfzLJisG8R\nERHpoUi6IcI1jsCVhQbgPCAXeAQ4AJjR2QZmNhg4HVgT3E5ERERCkwWMAV52zm2JxgvGIiykAa3A\npc65nQBmdhPwrJld65xr7GSb04GnYlCbiIhIsroM+E00XigWYWEjsGFvUAj6BDBgFIEBjx2tAfj1\nr3/NoYce2usFxrOSkhLmzp3rdxlxQcciQMchQMdhHx2LAB2HgE8++YTLL78cgp+l0RCLsPAmcKGZ\n9XXO7Q4um0jgasP6LrZpADj00EMpKCiIQYnxKycnJ+WPwV46FgE6DgE6DvvoWAToOHxJ1LrxI5ln\noZ+ZHW5mRwQXjQv+Pjq4fo6ZLWqzyW+ALcDjZnaomZ0I3Ac81kUXhIiIiMSRSO6GOAr4ACgDHPAA\nUA7cGVw/DBi9t7FzbhdwGjAQeB/4FfA8cEPEVYuIiEjMRDLPwt/wCBnOuemdLFtBYNCiiIiIJBjN\nohjniouL/S4hbuhYBOg4BOg47KNjEaDj0HvMOed3DV9iZgVAWVlZmQariIiIhKG8vJzCwkKAQudc\neTReU1cWRERExJPCgoiIiHhSWBARERFPCgsiIiLiSWFBREREPCksiIiIiCeFBREREfGksCAiIiKe\nYvHUyaRWNLWIiuqKLtfn5+WzbPGyGFYkIiISXQoLPVRRXcGKohVdN1BOEBGRBKduCBEREfGksCAi\nIiKeFBZERETEk8KCiIiIeFJYEBEREU8KCyIiIuJJt072UH5efrvbI6t2VlHfVM+YQWP2rRcREUlg\nCgs91HHCpSc+fILpz0/n3VvfZWDWQJ+qEhERiR51Q0TZ5NGTAXh3/bs+VyIiIhIdCgtRNuGACQzO\nHszb69/2uxQREZGoUFiIMjPj2FHH8ta6t/wuRUREJCoUFnrB5NGTeXfDu7S6Vr9LERER6TGFhV4w\nadQktjdu5+PNH/tdioiISI8pLPSCo0ceTZql8fY6jVsQEZHEp7DQC/r36c/hQw/nrfUatyAiIolP\nYaGXTBo1SVcWREQkKYQdFszsBDNbYmYbzKzVzM4JY9vjzKzJzMrD3W+imTR6Esu3LGfL7i1+lyIi\nItIjkVxZ6Ad8CFwLuFA3MrMcYBHwagT7TDiTRk0C4J317/hciYiISM+EHRaccy855253zj0PWBib\nLgCeAlLi03PcoHHk9cvT5EwiIpLwYjJmwcymA2OBO2Oxv3hgZoFxCwoLIiKS4Ho9LJjZBOAe4DLn\nUmuWokmjJvHu+ndpbm32uxQREZGI9epTJ80sjUDXwx3Ouc/3Lg51+5KSEnJyctotKy4upri4OHpF\n9qLJoyezq2kX/6r+F0cMO8LvckREJMmUlpZSWlrablldXV3U92POhTxG8csbm7UC5znnlnSxPgfY\nCjSzLySkBf/cDBQ55/7ayXYFQFlZWRkFBQUR1+e3+qZ6Bvy/ATw05SG+d/T3/C5HRERSQHl5OYWF\nhQCFzrmo3H3Y290Q24GvAkcAhwd/FgCfBv+c1M9xzs7M5ohhR2hyJhERSWhhd0OYWT9gPPuuFIwz\ns8OBWufcOjObA4xwzl3hApctPu6wfTXQ4Jz7pIe1J4TJoybz4soX/S5DREQkYpFcWTgK+AAoIzDP\nwgNAOfvudBgGjI5KdUlg0uhJfL71c6p3VftdioiISEQimWfhb865NOdceoefq4LrpzvnTvbY/k7n\nXOIORAiTJmcSEZFEp2dD9LIDcw5keP/hvLVO4xZERCQxKSz0MjNj8ujJmpxJREQSlsJCDEwaNYn3\nN7xPU0uT36WIiIiETWEhBiaNnkR9cz0fVX3kdykiIiJhU1iIgYLhBfRJ78Pb69QVISIiiUdhIQay\nMrIoGF6gcQsiIpKQFBZiZNKoSbojQkREEpLCQoxMHj2ZiroKNu7Y6HcpIiIiYVFYiJG9kzOpK0JE\nRBKNwkKMjBwwktEDRqsrQkREEo7CQgxNGj1JVxZERCThKCzE0ORRkymrLGNPyx6/SxEREQmZwkIM\nTRo9icaWRj7Y+IHfpYiIiIQsw+8CUkXR1CLWVK2BWjj7+bMZlD2o3fr8vHyWLV7mT3EiIiIeFBZi\npKK6gpWnrwRgc/C/dpQTREQkTqkbQkRERDwpLIiIiIgnhQURERHxpLAgIiIinhQWRERExJPCgoiI\niHjSrZNR5JyjpOQ25s69GzNrty4/L7/d7ZFb67eyeddmxg0aR0Z6RmC9iIhIHFJYiKKysjLmzfs5\n3/72BRQWFrZb13HCpW0N2xj54EiuPO5KZn9jdizLFBERCYu6IaJo3rxnaW5+gEceebbbtgOzBnLZ\nf1zGwrKFNLc2x6A6ERGRyCgs9NDtt88hL28iEyacQWnpCuBqli5dzvjxU8jLm8jtt8/pcttrj76W\nDTs28MLyF2JXsIiISJgUFnpo9uybmT37FnbsyKCxcTFgVFUtZufOTGbPvoXZs2/uctsjhh3BpFGT\nmP+P+bErWEREJEwKCz2UmZnJrFkzyMlpP6AxJ8eYNWsGmZmZntt/76jv8eqqV1mxZUVvlikiIhKx\nsMOCmZ1gZkvMbIOZtZrZOd20n2pmy8ys2szqzOwtMyuKvOT41NraTHb2o+TknAY8SmNjaOMQLjrs\nIgZnD2bBPxb0boEiIiIRiuTKQj/gQ+BawIXQ/kQCNw2eARQAfwFeMLPDI9h33CooGMu99xqfffY8\n/foZmZljQ9ouKyOLGUfO4PEPH2d30+5erlJERCR8YYcF59xLzrnbnXPPAxZC+xLn3E+cc2XOuc+d\nc7cBK4GzI6g3bj399DxmzZpBbm5fbr55BpWV89iyJbRtrznqGuoa6vjtv37bu0WKiIhEIOZjFiww\nW9H+QG2s9x0r114Lra2wIMSehXGDxnHGhDOY9/48nAvlYo2IiEjs+DHA8YcEujKe8WHfMZGXB1dc\nAT//OTQ0hLbNtUddS/nGct6vfL93ixMREQlTTMOCmV0KzAYucs7VxHLfsXbTTVBdDb/+dWjtp4yf\nwpiBY5j/vm6jFBGR+GI9uextZq3Aec65JSG0vQR4FLjQOfdSN20LgLITTzyRnJycduuKi4spLi6O\nuOZYmjoVPv0U/v1vSAshlt37xr3c8dc72HDTBgb3Hdz7BYqISEIrLS2ltLS03bK6ujr+/ve/AxQ6\n58qjsZ+YhAUzKyYQFC52zi0N4XULgLKysjIKCgoirs9vb74Jxx8PL7wAZ53VffvNuzYzau4o7jn5\nHn4w+Qe9X6CIiCSd8vLyvc8nilpYiGSehX5mdriZHRFcNC74++jg+jlmtqhN+0uBRcAPgPfNbGjw\nZ0A03kA8mzwZjj0WfvKT0NoP6TeEi75yEY/84xFaXWvvFiciIhKiSJ46eRSBuRJc8OeB4PJFwFXA\nMGB0m/ZXA+nAvOAPHdonLTO4+Wa48EJ4/304+mjv9kVTi/h0w6esq1vHgb86kH59+rVbn5+X/6Wn\nV4qIiPS2sMOCc+5veFyRcM5N7/D7SRHUlTTOOw8OOihwdeHpp73bVlRXsO7MdQBsYMOXGygniIiI\nD/RsiF6Wnh64M+J3v4PVq/2uRkREJHwKCzFw5ZUwaBD89Kd+VyIiIhI+hYUY6Ns3MKvjY49BbdLO\nWykiIslKYSFGrr8empth4UK/KxEREQmPwkKM7J0C+qGHoKHBceON/6XnQIiISEJQWIihm26CTZvg\nnnvKmDfv55SXhzdXRmur5l4QEZHYi2SeBYnQxIlwzjnw0EPP0tz8AI888iyPPlr4xfr8vPxOb49s\nbG5kbd1aagfW0tTSRGZ6ZgyrFhGRVKewECO33z6HBQueoE+fcdTVZQH/j6VLz2f8+Cls376amTOv\n9Jxw6c+r/8yUX0/hmqXX8Ng5jxF40reIiEjvUzdEjMyefTOzZ99Cc3MGsBgwqqoWs3NnJrNn38Ls\n2Td7bn/y2JP55bm/5PEPH+euv90Vk5pFRERAVxZiJjMzk1mzZvDww89TVbVveXa2MWvWjJBe4/Kv\nXc7aurXcdv1tPJLxCDlZOZ2207TQIiISTQoLMdba2kx29qMMGPA01dUXs2FDM1u3BiZtCsWPjv8R\n93EfVedUUUVV542UE0REJIrUDRFjBQVjufdeY9Wq57ntNgPGMm0aNDWFtr2Zkdcvr1drFBERaUtX\nFmLs6af3PXjzxz+ewcknz6CoCG68EebN89iwDQ1uFBGRWNKVBZ+ddBLMnx/4CTUsiIiIxJKuLMSB\nq6+Gjz+GG26ACROgqKhnr1fXUMep553Kus3rumyjQZAiIhIqhYU4cf/9sHw5TJsG77wDEyc6Skpu\nY+7cu8PudqjaWUXN8hpaLmnpupFygoiIhEhhIU5kZEBpKUyeDGedBQsWBKaE/va3L6CwsLD7F2hj\nzMAx1NbXsp3t3bYtmlpERXVF4BcH1ZtryBuSC8F8oisQIiKisBBHcnLghRfgmGPg29/ufEpoaD8t\ndENDA2vXruPAA0eTlZUVWD8in4rqCs+wsPchVhXVFawoWtFu3TbaPEdbOUFEJOVpgGMcuf32ORx7\n7ET22+8MNm1aAVxNaely8vKmMGjQRG64YQ6trbBs8TKWv7mc5W8u59QjLoGaBZx2ZPEXy/747DKq\nq2s997Vyy0pOW3QGtfVbY/PmREQkYSksxJHOpoTevXsxmzdnsm3bLTz00M307QtDhsxhv/0mMnDg\nGTz7bCBUPPXUcrKzp5CePpE+feawbWs3szztzuXV11qo2bnZs1lNjXfoEBGR5KewEEf2Tgmdk9N+\nQOOECcby5TN48cVM7rsPpk27mYMOuoWdOzPYsSMQKhoaFgOZnHzyLTz88M2MHOU9KHJ07gE8fdYy\nctPHe7YbPHgQRVOLmHjcxMDP5IkMmjCYiZMnfrGsaGoPb98QEZG4pjELcWjvlNADBz7Ntm0X41wz\nBx8MBx+8t0UmMIOJE59nRZvhBgceaLzySuA5Ez97ynsf2X0Dd17M/lkaNR7t1myuZmf9HjaeV9Fu\nucY1iIikDoWFOFRQMJbvf9+YMeN5HnuslDfeGNtpu46horW1+Yt1ewdBVlZWkZ2dxcCcAWyr2059\nfQMjRgwNrA9Bk+1k4+66qLwvERFJTAoLcajtlNCzZs3o8qmUXqEiWrc7Thg8jnXVm2hgR5dtNK5B\nRCS5KSwksFBDhZe2t2F2un5YPhis9AgLgweH+MhMERFJSAoLKS6UKxATj5vouV4PthIRSW66G0J6\nbEf9br9LEBGRXqSwID22cdd6/s9zj/hdhoiI9JKwuyHM7ATgh0AhMBw4zzm3pJttvgk8ABwGrAXu\nds4tCrta8YXXuIbmJse61gbu/uhaFt59N4My+wa6JfScCRGRpBHJmIV+wIfAY8AfumtsZmOApcB8\n4FLgVOBRM6t0zr0Swf4lxrr7gN+5EyZf/yj/3H41NZe2X6f5GEREEl/YYcE59xLwEoCFNrLte8Aq\n59wtwd+Xm9nxQAmgsJAE+veH8ke/y6BD72In6/wuR0REoiwWYxaOBV7tsOxlYFIM9i0xkpEBfdLq\nPdtoPgYRkcQUi7AwDKjqsKwKGGBm+8Vg/xIj3c23oPkYREQSk+6GkKjprldK8zGIiCSmWEzKtAkY\n2mHZUGC7c67Ra8OSkhJycnLaLSsuLqa4uDi6FYqIiCSg0tJSSktL2y2rq4v+83xiERbeBs7osKwo\nuNzT3LlzKSgo6JWiREREEl1nX6DLy8spLCyM6n4imWehHzCeL+6eZ5yZHQ7UOufWmdkcYIRz7org\n+gXAdWZ2L/BL4BTgQuBbPa5e4kpn8zHs3NlK5a4K0jPSGH3IaH8KExGRHonkysJRwF8AF/x5ILh8\nEXAVgQGNX3wqOOfWmNmZwFzg+8B6YIZzruMdEpLgupqPofimcn7b/1hGHXJkjCsSEZFoiGSehb/h\nMTDSOTe9k2V/JzDjo6Sgx+8p4K+XzmFR+s1c9EkRZx56mt8liYhIGHQ3hPS6rCx46f+WYKuKmPbb\n77B512a/SxIRkTAoLEhMHP61NO48YhG761s469HpOOf8LklEREIUi7shRAC47YZh3H/wSN5Lf5G8\nJ4ZxQN+BX2qjh02JiMQfhQWJmbQ0GDJ4Nzu+BTVUU0P1lxspJ4iIxB11Q0hMZSieiogkHIUFERER\n8aSwIDHV3ZMn9WRKEZH4o7AgMaUnU4qIJB6FBYkpPZlSRCTxKCyIiIiIJ41Nl5ja+7CpysoqsrOz\nyMkZwKrKjbRm7WTsoHHkD8/3u0QREelAYUFiqrMJl5770w6m/nUMXx17BktmPuxDVSIi4kXdEOK7\nc6fsz+j1N7G08lE2bK/0uxwREelAYUF8ZwYPXnI9bk823//t/X6XIyIiHSgsSFy44KwcRqz7Ps+t\nX8imHVV+lyMiIm0oLEhcMIOfXHQDrc0Z3PjMA36XIyIibSgsSNy45NwDGFpxPc+umU/N7hq/yxER\nkSCFBYkbZnD/+TfR2golT//U73JERCRIYUHiyuXn5zKk4nuUrnqIrfVb/S5HRERQWJA4YwZzzv4B\nLa6Jm599yO9yREQEhQWJQ9MvGsbgNdfw5Mqfsr1xu9/liIikPIUFiTtpaXD3t35IM7u55Xea0VFE\nxG+a7lni0u+enU7aR9kszJjNaw8/Tk1NLXlDciH4UMr8vPxOp44WEZHoU1iQuLS2uoLWK+oA+IzP\nANhG7b4GygkiIjGjbggRERHxpLAgcammprZH60VEJHoUFiQuDR48qEfrRUQkehQWJC6ZWY/Wi4hI\n9EQUFszsOjNbbWb1ZvaOmR3dTfvLzOxDM9tlZpVm9piZHRBZySIiIhJLYYcFM7sYeAC4AzgS+Ah4\n2cxyu2h/HLAI+B/gK8CFwNeBX0RYs4iIiMRQJLdOlgALnXNPApjZTOBM4Crgvk7aHwusds7NC/5e\nYWYLgVsi2LekiPy8/Ha3R9Zsa6C2eS3D+o1kQHa/wHoREYmJsMKCmWUChcA9e5c555yZvQpM6mKz\nt4G7zewM59yfzGwocBHwYoQ1SwroOOHSjh2OnP86jCFDjuB/b/+NT1WJiKSmcLshcoF0oKrD8ipg\nWGcbOOfeAi4HnjazPcBGYCtwfZj7lhS2//5GQdqV/KtpMVvrt/ldjohISun1GRzN7CvAz4D/S+DC\n8nDgJ8BC4Lte25aUlJCTk9NuWXFxMcXFxb1Sq8S3H515ORe++SPuef4Z7r/kP/0uR0TEd6WlpZSW\nlrZbVldXF/X9mHMu9MaBbojdwAXOuSVtlj8B5DjnpnayzZNAlnNuWptlxwGvA8Odcx2vUmBmBUBZ\nWVkZBQUFYbwdSWatrdD/mjPYP3c7VXPe9LscEZG4VF5eTmFhIUChc648Gq8ZVjeEc64JKANO2bvM\nAje8nwK81cVmfYHmDstaAccXjwUS6V5aGpw16kqqs97ig7Ur/C5HRCRlRDLPwoPA1Wb2HTM7BFhA\nIBA8AWBmc8xsUZv2LwAXmNlMMxsbvKrwM+Bd59ymnpUvqebHl58LDTnc/odF3TcWEZGoCHvMgnPu\nmeCcCncBQ4EPgdOdc5uDTYYBo9u0X2Rm/YHrCIxV2Aa8Bvx/PaxdUtDEg7IYWVvMKzxJS+tdpKel\n+12SiEjSi2gGR+fcfOfcGOdctnNuknPuH23WTXfOndyh/Tzn3H845/o750Y5565wzm3safGSmr57\n9BU0Zq3nN2//xe9SRERSgp4NIQnnh8XHkFY7kftfecLvUkREUoLCgiScfv2MwvQr+VfLH9i6O/q3\nCImISHsKC5KQ/uvMb+PSGvnv5571uxQRkaSnsCAJ6dyTR9K36lR+/b+6K0JEpLcpLEhCMoOzR11J\ndfYbfLD2M7/LERFJagoLkrB+fPl50DCA23+nqwsiIr1JYUES1oQx2Yyqu4Rl1Ytoda1+lyMikrQU\nFiShXX30lezJXsevXtecCyIivaXXnzop0luKphaxZlMF1GTyn7+Zyo2tmeQNyf3iiSP5efksW7zM\n3yJFRJKAwoIkrIrqClZOCTxQag9N7AG2UbuvgXKCiEhUqBtCREREPCksSMKqqant0XoREQmNwoIk\nrMGDB/VovYiIhEZhQRKWmfVovYiIhEZhQURERDwpLIiIiIgn3TopCSs/L7/d7ZG7djezYfcqBvUZ\nxpABAwLrRUSkxxQWJGF1nHDJOciadQyDBoxi+T2/96kqEZHko24ISRpmcFS/qaxOe4ldjfV+lyMi\nkjQUFiSpzPzmVFzmbuYv0/SNIiLRorAgSaX4tImkb/kKT77/B79LERFJGgoLklQyMuCwjKl83PwC\nTS1NfpcjIpIUFBYk6Vx5zFRa99vKr17/u9+liIgkBYUFSTozzynA6g7kF6+rK0JEJBoUFiTpZGcb\nBzVNpXz3c7S6Vr/LERFJeAoLkpQuOXwqTVmVLC1/3+9SREQSnsKCJKWSC46HXbn8/NXFfpciIpLw\nFBYkKR0wKJ2RO8/lzdo/4JzzuxwRkYQWUVgws+vMbLWZ1ZvZO2Z2dDft+5jZ3Wa2xswazGyVmV0Z\nUcUiITpv4lTq+67krZUf+12KiEhCCzssmNnFwAPAHcCRwEfAy2aW67HZs8BJwHTgYKAYWB52tSJh\nuPmCU6CxPw/+SV0RIiI9EcmVhRJgoXPuSefcp8BMYDdwVWeNzWwKcALwLefcX5xza51z7zrn3o64\napEQjBmVxeDaM3mtUrdQioj0RFhhwcwygULgtb3LXKBD+FVgUhebnQ38A7jVzNab2XIzu9/MsiKs\nWSRkU8ZMpa7vB3xcucbvUkREEla4VxZygXSgqsPyKmBYF9uMI3Bl4TDgPOAG4EJgXpj7FgnbD849\nA5r78JOlz/ldiohIwsqIwT7SgFbgUufcTgAzuwl41syudc41drVhSUkJOTk57ZYVFxdTXFzcm/VK\nEjny0AH033wqLzYtBm70uxwRkagqLS2ltLS03bK6urqo7yfcsFADtABDOywfCmzqYpuNwIa9QSHo\nE8CAUcDnXe1s7ty5FBQUhFmiSHvfHHo+S+1qNmyrZuTAPL/LERGJms6+QJeXl1NYWBjV/YTVDeGc\nawLKgFP2LjMzC/7+VhebvQmMMLO+bZZNJHC1YX1Y1YpEoOSMc8AZDyxd4ncpIiIJKZK7IR4Erjaz\n75jZIcACoC/wBICZzTGzRW3a/wbYAjxuZoea2YnAfcBjXl0QItFy0jFD2K/6eH7/sW6hFBGJRNhh\nwTn3DHAzcBfwAfA14HTn3OZgk2HA6DbtdwGnAQOB94FfAc8TGOgo0uvM4NgB57M241W21W/3uxwR\nkYQT0QBH59x8YH4X66Z3smwFcHok+xLpqaKpRaxY/RnU72H8cxNpqd9D3pDcwKgZID8vn2WLl/lb\npIhIHIvF3RAivqqormDj1NUAbAmOw91G7b4GygkiIp70ICkRERHxpLAgSa+mprZH60VEUp3CgiS9\nwYMH9Wi9iEiqU1iQpBeYCiTy9SIiqU5hQURERDwpLIiIiIgn3TopSS8/L7/d7ZGfbdoI6Y2MHzJm\n33oREemSwoIkvY4TLl3y46d5uvUSln7vFSbkHehTVSIiiUPdEJJyZn1rCrRk8LOXXvC7FBGRhKCw\nIClnckEO+1V9gxeWKyyIiIRCYUFSjhkcPeBs1qb/he0NO/wuR0Qk7iksSEq66rizIX0Pj/1ND4YQ\nEemOwoKkpOIp47DNh/Gr95b4XYqISNxTWJCUlJUFE1rP4Z8NL9LS2uJ3OSIicU1hQVLWBf9xNs19\ntvDHf73tdykiInFNYUFS1vXnfR125rHgNd0VISLiRWFBUtaI4ekM2XYmr1dr3IKIiBeFBUlppx14\nDjuyPuXjTZ/5XYqISNxSWJCUNutbp0Hzfjz0sroiRES6orAgKe2Ygn5kbTyFpSvVFSEi0hWFBUlp\nZvD1nLPZkP46tbu3+l2OiEhcUliQlPfdE8+CtBYe/duf/C5FRCQuKSxIyrvo9FGkbSrgqX9o3IKI\nSGcUFiTlZWXBwe4c/t34J5pamvwuR0Qk7igsiADTjjiblsw6lv7zdb9LERGJOwoLIsDMc4+E7SNZ\n8BfdFSEi0pHCgggwfLiRt/Vs3qxZgnPO73JEROJKRGHBzK4zs9VmVm9m75jZ0SFud5yZNZlZeST7\nFelNp489h119VvNR5cd+lyIiElfCDgtmdjHwAHAHcCTwEfCymeV2s10OsAh4NYI6RXrd9WeeBHv6\n8fAruitCRKStjAi2KQEWOueeBDCzmcCZwFXAfR7bLQCeAlqBcyPYr0ivKZpaREVVBVYFT6Tfyev/\n83i79fl5+SxbvMyn6kRE/BVWWDCzTKAQuGfvMuecM7NXgUke200HxgKXAbMjK1Wk91RUV7Di9BUA\ntAArWNG+gXKCiKSwcK8s5ALpQFWH5VXAxM42MLMJBMLF8c65VjMLu0gRERHxTyTdECEzszQCXQ93\nOOc+37s41O1LSkrIyclpt6y4uJji4uLoFSkiIpKgSktLKS0tbbesrq4u6vuxcG4TC3ZD7AYucM4t\nabP8CSDHOTe1Q/scYCvQzL6QkBb8czNQ5Jz7ayf7KQDKysrKKCgoCOf9iERk8MQh1F5a0+X6A36T\ny5blm2MXlb8hAAATmElEQVRYkYhIZMrLyyksLAQodM5F5e7DsO6GcM41AWXAKXuXWaBf4RTgrU42\n2Q58FTgCODz4swD4NPjndyOqWiTKBg8e1KP1IiLJLJJuiAeBJ8ysDHiPwN0RfYEnAMxsDjDCOXeF\nC1y2aHfTuplVAw3OuU96UrhINHU3lkZjbUQklYUdFpxzzwTnVLgLGAp8CJzunNt7jXYYMDp6JYqI\niIifIhrg6JybD8zvYt30bra9E7gzkv2K9Jb8vHxYBpWVVWRnZ+Esi5rmCvZr6kv+8FGB9SIiKapX\n74YQSRQdJ1yqq4NBtx7FIaPH8uFtz/pUlYhIfNCDpEQ6kZMDh7ZO45+NL7Jzz06/yxER8ZXCgkgX\nrj7uIlrT63nirRf9LkVExFcKCyJdmHH+WKzyaBa+8YzfpYiI+EphQaQL++8Ph9lF/HvPH9UVISIp\nTWFBxMM1x1+ES2/gsdeX+l2KiIhvFBZEPEyfOoa0yq/zP2+pK0JEUpfCgoiHfv3gP9Km8XHTH9nR\nuMPvckREfKGwINKNmSdeiEtv5Bd/e8HvUkREfKGwINKNK87LJ63yGB59W10RIpKaFBZEupGdDYdn\nTGN5y0tsb9zudzkiIjGnsCASgmu/GeiKWPAXdUWISOpRWBAJweVnH0h65SR++a66IkQk9SgsiIQg\nKwuO7DONFS0vsa2+zu9yRERiSmFBJETXnXQhLn0Pj/x5id+liIjElMKCSIguPWsU6ZWTefw9dUWI\nSGpRWBAJUZ8+ULjfND5zy9hav83vckREYkZhQSQMs04NdEU8/Iq6IkQkdSgsiITh4jNGkrHheBb9\nQ10RIpI6FBZEwpCZCUf3ncbnLKN291a/yxERiYkMvwsQSSRFU4tYs+pzaGhiwnOH0lrfRN6QXLDA\n+vy8fJYtXuZvkSIiUaawIBKGiuoKNp6/CoBaqgDYRu2+BsoJIpKE1A0hIiIinhQWRMJQU1Pbo/Ui\nIolIYUEkDIMHD+rRehGRRKSwIBIGM+vRehGRRKSwICIiIp4iCgtmdp2ZrTazejN7x8yO9mg71cyW\nmVm1mdWZ2VtmVhR5ySIiIhJLYd86aWYXAw8A/wm8B5QAL5vZwc65mk42OZHADWU/ArYBVwEvmNnX\nnXMfRVy5iA/y8/Lb3R7Z0gKf16whIz2TcbkjA+tFRJJMJPMslAALnXNPApjZTOBMAiHgvo6NnXMl\nHRbdZmbnAmcDCguSUDqbcOns237F0j7f4ZeX/Z7jxn/Vh6pERHpXWN0QZpYJFAKv7V3mnHPAq8Ck\nEF/DgP0B3WMmSeGR6y6ButFc99SXsrKISFIId8xCLpAOwanr9qkChoX4Gj8E+gF6Eo8khVEjMjk+\n/SY+ainlk8q1fpcjIhJ1Mb0bwswuBWYDF3UxvkEkIS285ruwZ3+uWfSg36WIiERduGMWaoAWYGiH\n5UOBTV4bmtklwC+AC51zfwllZyUlJeTk5LRbVlxcTHFxccgFi8TCV8b35/DG63lj1wNUbpvNiIGD\n/S5JRFJAaWkppaWl7ZbV1dVFfT8WGHIQxgZm7wDvOuduCP5uwFrgIefc/V1sUww8ClzsnFsawj4K\ngLKysjIKCgrCqk/EL299uJnjfpfPOQf8iOdvmu13OSKSosrLyyksLAQodM6VR+M1I+mGeBC42sy+\nY2aHAAuAvsATAGY2x8wW7W0c7HpYBPwAeN/MhgZ/BvS4epE4MvmIIRy0/SqWbn6I7fW7/S5HRCRq\nwg4LzrlngJuBu4APgK8BpzvnNgebDANGt9nkagKDIucBlW1+fhp52SLxae60H9DaZys3PPG436WI\niERNJPMs4JybD8zvYt30Dr+fFMk+RBLR2cePZfjT03iq/icsbL6GPhkR/RUTkSgqmlpERXVFl+v3\nTqbWXZvO5llJFfqXTCTK/vtbP2TGewX811PP8pMrNBhXpLeEEgKWLV5GRXUFK4pWdP1CwQwQSptU\npbAgEmVXnXEkP/hTEY9sv5f7Wi8hLU1PohQJVyhBINQQID2nsCASZUVTi8j47BO27VnP0FdG01xf\nT96QXAhmhlS/nCkSimgEAeccG3dspLG50bNdze4awr0zMNUoLIhEWUV1BTUXrgeghg0AbGs7u7ly\ngqSwULsOurOneQ9NrU2ebVZuWcmIB0cEHmHooa5+O62trZ5tampS+wkFCgsiIhIzoVwxcM7R0tri\n+Tprtq3pdl+2qz99ny9lV93lQNcTFbXUDAQGASu7bDN48KBu95fMFBZEoqy7byCp/g1FklO7KwYO\nqjfXRNT9trZuLUPuH8KW2i3eDbfnQMsBwOoumwzsO5xbLz6LBxflUe0RFsaNG4SZ8bnH7gLzD6Yu\nhQWRKBs8eBC1dP3ok1T/hiLJqbMrBu26316Gz2s/Z0fjDs/XadydyZ4/3Qg75gMbu2w3YkgeGZmG\n16Pbhgwxbr0VfrnEqPZol5GZ2kEgFAoLIlHW3TeQVP+GIoknlKsG3VlZ+xnjfz4evLMC++3Zj7kX\n/x/uXvAr1nuEhf77h/73KD8vH5bB6lXraW7en7T0HbS27E9Gxg7Gjhu1r/5lUFlZRXZ2FgNzBrCt\nbjv19Q2MGDE0pPeYzBQWRETEU3dXDdzLjqotXt/dwe3Kwp76Pem7rqOZVV22G33gIGbOhLm/6r6u\nvSGguw/4vd0fF198HccfX8CMGcU89lgpb7xRztNPz+t+R6KwIBJr9Y17/C5BBIjenQkrt3wGeN96\neOCQkXzyxhSOPC0Dj+GNX1x5CyUIhHsLcttgMGvWDGbNmhHW9qlMYUEkyvb+I9fRnqYW1mxdy/r9\nqvh8UzUHDcuLfXEibYQ6l8Hmmm4GG+7Opk9aLns8RhBkZRt9+4Zem+YiiS8KCyJR5vWP3Cvvr+H0\nZyZR+MBZVNz5F3L69othZSLhqdpew0H3TGYr3mFh/KiRpKWZ5xWDvULtOpD4orAgEkOnHT2G+Rv+\nyPfeO5HDf3wJK/97MZnp+mso0ReNQYl1jduo+3AUGY15NHvcTxDOlOa6YpCY9K+USIzNPO9IPl3/\nO372iykMWJzLgQcMBSK/L12kM90NSmz+UzMbuxmUuH/zIDY89gxHTZnICs+bD7vufmu3XhKWwoKI\nD356/en8z0ND2X1ZFSvaTBajaaElVlZtXQVp3m2GDTuA/fcP7fUUbJObwoKIT0bl5rCCKr/LkAQT\n6kyJ3c0Umtm4P7k5eWz0mLew450JXdFVg+SnsCDik5otmhZawtdd90Lt73dx9B2zqHVdT28MMGbE\nMMzMY9qjfXTVQBQWRHyiaaGlo2gMSqxp3EDtjqVkur40eTwPwcx0Z4KETGFBxCfdTfu8c88uTpt6\nGmurg/eu9+DhPJIYun2+wjKo3ux9G2P/loFsvW8Vh717SLvxMJ3RuSOhUlgQiVMbd1RS80ktTcUN\n7ZZrEGTq+qxqPa1puz3bDM8bQkaGnj8i0aWwIBKnsptHUu/RTSGJI1qDEh176J+ey06P80KDEqU3\nKCyI+KTbf8y/ls+7//sB22nsso0GQSaG7roXWl+ChX/4hN3ej1dgfN5YzIwVIYRIdTFINCksiPgk\nlH/MD558MNs9PhgGDOwPhP7NVaIvGoMSP6tZxcx/fuWLbbqix5uLXxQWROJYdx8Oa7atYcwdJ7Nr\n1SfUnL+p3TqNbYiNUAYldncFKKM1g/knvsB9L3+fz1jp2VbdC+IHhQWRBJbdMoz1la20NGzqvrH4\noqK6lsZW70GJY4eO5uqTpvCTEK4c6AqR+EFhQSSBjR46gH/O+ysHHHIAu9jaZbuamlp1VYQpWoMS\nG62GjPS+NHu00aBEiXcKCyJxLJQPjz59YMSwXFZ6hIWd7OHD1f9m89TKdsvVVQHOOUpKbmPu3Lvb\ndft0173Q+ILj1offZleL9+uPHzyetLQ0VoTwAGeFNYlXCgtxrrS0lOLiYr/LiAupeCw6+/Do7Dh0\nN7ZhT0s9m+u3h7TPRLkC0ZPzoe17bGhoYO3adSx+/SmysrKA0L7BV+xYyX1bJkO697FPS+vmaU1R\nkIp/Nzqj49B7IgoLZnYdcDMwDPgImOWce9+j/TeBB4DDgLXA3c65RZHsO9Xo5N9HxyIgkuMwYchY\nNlRtZrfH1Yfqms0450IasNcuUHQiFoGis+MQatDp7D2uZe0Xf979XCubt3k/knm/ln68fOlfufrl\nS1np86BE/d0I0HHoPWGHBTO7mMAH/38C7wElwMtmdrBz7kv3eJnZGGApMB+4FDgVeNTMKp1zr0Re\nuoiEytJgZDddFdvYSt/bDqS1xnuKYOj8En07wQ/GUG8rDOUDvuNrrfpoNRMnT2zXLpSgA93fnbB+\n92fQx7MJB+aN4BsTjgrpdsZ4uBIj0hORXFkoARY6554EMLOZwJnAVcB9nbT/HrDKOXdL8PflZnZ8\n8HUUFkSiIJRvrl5XAgD2ZwS5Wy5hdfNPPdtV12wmb/CQkOoK9cM7lDZfeq0tsOL0FV9q56WqdieT\nf3g/dc17PNsdsN9Q+mdlsZauj5kGJUoqCSssmFkmUAjcs3eZc86Z2avApC42OxZ4tcOyl4G54exb\nRLoWyjfXicdN9Fw/PLc/yxfez8GTnve8rL6NrWwL4epDrG3ftYdt272vGNS1VvJu1p2Q7nVvAuTu\nnxPyfnXVQFJBuFcWcoF0oKrD8iqgq3+JhnXRfoCZ7eec62wu2yyATz75JMzykk9dXR3l5eV+lxEX\ndCwCIj0ODTsboNJ7fXl5OY27Gz3b9duTS33jTlorG7pss2LdCrKuOJDGtZWer7W6Yh2GebZZtWYt\nBd+dyWdrKtq3a6Dd75u2rgn8weO1RmSOYMl5Szh/yfmsrVzbZbuGnQ3dvtbe4xUP9HcjQMchoM1n\nZ1a0XtOc62Yy8raNzYYDG4BJzrl32yy/FzjROfelqwtmthz4pXPu3jbLziAwjqFvZ2HBzC4Fngrn\njYiIiEg7lznnfhONFwr3ykIN0AIM7bB8KNDVFHKbumi/vYurChDoprgMWEPg+4OIiIiEJgsYQ+Cz\nNCrCCgvOuSYzKwNOAZYAWGCUzynAQ11s9jZwRodlRcHlXe1nCxCVNCQiIpKC3ormi0UyW8iDwNVm\n9h0zOwRYAPQFngAwszlm1nYOhQXAODO718wmmtm1wIXB1xEREZE4F/atk865Z8wsF7iLQHfCh8Dp\nzrnNwSbDgNFt2q8xszMJ3P3wfWA9MMM51/EOCREREYlDYQ1wFBERkdTT+5OWi4iISEJTWBARERFP\nvoQFM7vOzFabWb2ZvWNmR3fT/ptmVmZmDWa2wsyuiFWtvS2cY2Fm3zCz1g4/LWaWF8uao83MTjCz\nJWa2Ifiezglhm6Q7J8I9Dkl8PvzIzN4zs+1mVmVmi83s4BC2S8ZzIuxjkYznhZnNNLOPzKwu+POW\nmU3pZpukOx8g/GMRrfMh5mGhzYOo7gCOJPDUypeDgyY7az+GwAROrwGHAz8j8CCq02JRb28K91gE\nOWACgYGkw4Dhzjnvx+PFv34EBspeS+D9eUricyKs4xCUjOfDCcDPgWMIPHguE1hmZtldbZDE50TY\nxyIo2c6LdcCtQAGBRw78GXjezA7trHESnw8Q5rEI6vn54JyL6Q/wDvCzNr8bgTskbumi/b3A/3ZY\nVgr8Mda1x8Gx+AaBSbEG+F17Lx6TVuCcbtok7TkR5nFI+vMh+D5zg8fj+FQ+J8I4FqlyXmwBpqfy\n+RDisYjK+RDTKwttHkT12t5lLvBuInkQVVftE0KExwICgeJDM6s0s2VmNrl3K41LSXlORCgVzoeB\nBL4ZeT0lKlXOiVCOBSTxeWFmaWZ2CYH5fbqa3C8lzocQjwVE4XyIdTeE14OohnWxjeeDqKJbXkxF\nciw2AtcAFwDnE7gc9VczO6K3ioxTyXpOhCvpzwczM+CnwBvOuY89mib9ORHGsUjK88LMvmpmO4BG\nYD4w1Tn3aRfNk/p8CPNYROV8CHtSJvGPc24FsKLNonfM7CCgBEiKwTsSuhQ5H+YDXwGO87uQOBDS\nsUji8+JTAuMPcgjMAvykmZ3o8SGZzEI+FtE6H2J9ZSFWD6JKBJEci868B4yPVlEJIlnPiWhImvPB\nzB4GvgV80zm3sZvmSX1OhHksOpPw54Vzrtk5t8o594Fz7jYCA8Jv6KJ5Up8PYR6LzoR9PsQ0LDjn\nmoC9D6IC2j2IqquHXrzdtn2Q54OoEkGEx6IzRxC4zJRKkvKciJKkOB+CH47nAic559aGsEnSnhMR\nHIvOJMV50UEa0FWXQtKeD13wOhadCf988GHU5jRgN/Ad4BBgIYGRnEOC6+cAi9q0HwPsIDC6dSKB\n28r2AKf6PQLVh2NxA3AOcBBwGIH+yyYC3zZ8fz89OA79CFxSO4LASO8bg7+PTqVzIoLjkKznw3xg\nK4HbBoe2+clq0+aeFDknIjkWSXdeBN/jCUA+8NXg34Vm4OTg+pT4NyLCYxGV88GvN3stsAaoJ5D0\njmqz7nHgzx3an0jgW3g9sBL4tt//w/w4FsAPg+9/F7CZwJ0UJ/r9HqJwDL4R/HBs6fDzy1Q6J8I9\nDkl8PnR2DFqA77RpkyrnRNjHIhnPC+BRYFXw/+0mYNneD8dUOh8iORbROh/0ICkRERHxpGdDiIiI\niCeFBREREfGksCAiIiKeFBZERETEk8KCiIiIeFJYEBEREU8KCyIiIuJJYUFEREQ8KSyIiIiIJ4UF\nERER8aSwICIiIp7+f/M4H0vHR3oeAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f8fa4992a10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lattice = quippy.farray(np.eye(3)*20.0)\n",
    "at = quippy.Atoms(n=3,lattice=lattice)\n",
    "at.set_atoms([14 for i in range(at.n)])\n",
    "at.pos[:,1] = quippy.farray([0.0,0.0,0.0])\n",
    "at.pos[:,2] = quippy.farray([2.2,0.0,0.0])\n",
    "\n",
    "r = 3.0\n",
    "\n",
    "GAP3bEnergy = []\n",
    "SW3bEnergy = []\n",
    "angles = np.linspace(20.0,180.0,40) / 180.0 * np.pi\n",
    "for a in angles:\n",
    "    at.pos[:,3] = r*quippy.farray([np.cos(a),np.sin(a),0.0])\n",
    "    \n",
    "    pot.calc(at,energy=True,args_str=\"energy_per_coordinate\")\n",
    "    GAP3bEnergy.append(at.energy_per_coordinate[2])\n",
    "    potSW3.calc(at,energy=True)\n",
    "    SW3bEnergy.append(at.energy)\n",
    "    \n",
    "plt.plot(angles,GAP3bEnergy,\"*-\")\n",
    "plt.plot(angles,SW3bEnergy,\"s-\")"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [default]",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
