lionello/Coul_Buck.ipynb

## Coul_Buck.ipynb
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ao=1e-010\n",
    "eV=1.602176634e-019\n",
    "ke=8.987551787e+009 \n",
    "\n",
    "\n",
    "#Define a class \"Atom\" to carry the parameters of the atom that composes the lattice\n",
    "class Atom:\n",
    "    def __init__(self, name, sigma, latt_const, epsilon, charge, A, B, C):\n",
    "        #inputs\n",
    "        self.name = name #atom name\n",
    "        self.sigma = sigma #sigma of atom (distance where potential = 0) [for LJ]\n",
    "        self.latt_const = latt_const #lattice constant\n",
    "        self.epsilon = epsilon #depth of potential well for LJ\n",
    "        self.charge = charge #charge of the atom [coulumb]\n",
    "        self.net_charge =0 #sum of charge in the neighborhood\n",
    "        self.A = A #constant for Buckingham potential\n",
    "        self.B = B #constant for Buckingham potential\n",
    "        self.C = C #constant for Buckingham potential\n",
    "        self.pos = (0,0,0)\n",
    "        self.potential = 0\n",
    "        self.force_vec = (0,0,0)\n",
    "        self.dF_m=[[0,0,0],[0,0,0],[0,0,0]]\n",
    "        self.g = [0,0,1] #initialization of CG algorithm\n",
    "        self.h = [0,0,0] #initializatino of CG algorithm\n",
    "\n",
    "    def clone(self):\n",
    "        return Atom(self.name, self.sigma, self.latt_const, self.epsilon, self.charge, self.A, self.B, self.C)\n",
    "\n",
    "    def cloneAt(self, pos:tuple):\n",
    "        newp = self.clone()\n",
    "        newp.pos = pos\n",
    "        return newp\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "Na = Atom(\"Na\", 0 ,5.6402*Ao ,0 ,+0.988*eV ,7895.4*eV ,0.1709*Ao ,29.06*eV*Ao**6)\n",
    "Cl = Atom(\"Cl\", 0 ,5.6402*Ao ,0 ,-0.988*eV ,1227.2*eV ,0.3214*Ao ,29.06*eV*Ao**6)\n",
    "NaCl = Atom(\"NaCl\", 0 ,5.6402*Ao ,0 ,0 ,2314.7*eV ,0.2903*Ao ,0 )\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "def Coul_Buck(r:float, vec_r:tuple, atom1:Atom, atom2:Atom, basis:Atom): #input distance, and the class \"Atom\" containing parameters of the atoms\n",
    "    A = atom1.A if atom1.name==atom2.name else basis.A #if the 2 atoms are the same, pick the atom's parameter\n",
    "    B = atom1.B if atom1.name==atom2.name else basis.B #if the 2 atoms are different, pick the basis parameter\n",
    "    C = atom1.C if atom1.name==atom2.name else basis.C\n",
    "    CB = ke*atom1.charge*atom2.charge/r + A*math.exp(-(r/B)) - C/(r**6)\n",
    "    F_CB = ke*atom1.charge*atom2.charge/(r**2) + (A/B)*math.exp(-(r/B)) - 6*C/(r**7)\n",
    "    dF_CB = -2*ke*atom1.charge*atom2.charge/(r**3) - (A/(B**2))*math.exp(-(r/B)) + 42*C/(r**8)\n",
    "    F_vec=[0,0,0]\n",
    "    dF_vec=[[0,0,0],[0,0,0],[0,0,0]]\n",
    "    for i in range(len(vec_r)):\n",
    "        F_vec[i]=F_CB*(vec_r[i]/r)\n",
    "    for j in range(len(vec_r)):\n",
    "        for i in range(len(vec_r)):\n",
    "            dF_vec[j][i]=dF_CB*(vec_r[j]*vec_r[i]/(r**2))\n",
    "    assert round(F_CB)==round(math.sqrt(F_vec[0]**2+F_vec[1]**2+F_vec[2]**2))\n",
    "    return (CB,F_vec,dF_vec)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-2.2520399089347627e-28,\n",
       " [-2.2520399089347627e-28, -4.504079817869525e-28, -4.504079817869525e-28],\n",
       " [[4.504079817869525e-28, 9.00815963573905e-28, 9.00815963573905e-28],\n",
       "  [9.00815963573905e-28, 1.80163192714781e-27, 1.80163192714781e-27],\n",
       "  [9.00815963573905e-28, 1.80163192714781e-27, 1.80163192714781e-27]])"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Coul_Buck(1,(1,2,2),Na,Cl,NaCl)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x11fb5ab70>"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAEQCAYAAACQip4+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGAZJREFUeJzt3XtwXGd5x/Hfo9XKkmI7VrDUpNjGKCQIzJALKpfQcZUAxQnUFFqm6QADJYxNJ9x6odChE1r/0dLpDAWmwNhNU8oEEoYQOi5Dwi1RBQUH5JCLjZULCmDHEAmj2HIkW9Lu0z/OWXu9XmnPWnt2z9n9fmYUrc45Onp0vPntq3ff877m7gIApEdbowsAAFSH4AaAlCG4ASBlCG4ASBmCGwBShuAGgJSJLbjN7BYzmzCzfTU6391m9rSZfa1k+6vM7H4ze8DMvmdmz6vFzwOApIqzxf05SVtqeL5/kfS2Mts/K+kt7n65pC9K+rsa/kwASJzYgtvdRyT9pnibmV0ctpz3mtl3zWygivN9R9J0uV2SVoePz5d0+FxrBoA0aK/zz9sl6d3u/piZvUzSZyRds8xzvkvS181sVtIxSS9f5vkAINHqFtxmtlLSVZK+bGaFzSvCfW+StKPMtz3p7q+tcOq/kHSdu99nZh+U9HEFYQ4ATameLe42SU+HfdFncPc7Jd1Z7QnNrFfSZe5+X7jpS5LuXlaVAJBwdRsO6O7HJD1hZm+WJAtctszTTkk638wuDb9+jaQDyzwnACSaxTU7oJndJmlI0lpJT0n6qKR7FIwCuUhSVtLt7l6ui6Tc+b4raUDSSklHJN3g7t8wszcq6GbJKwjyd7r7eG1/GwBIjtiCGwAQD+6cBICUieXNybVr1/rGjRvjODUANKW9e/f+2t17oxwbS3Bv3LhRo6OjcZwaAJqSmf086rF0lQBAyhDcAJAyBDcApAzBDQApQ3ADQMrUe3bARQ2PTWjnyLgOTs1ofU+3tm/u19BAX6PLAoDESUSLe3hsQjft3q+J6RNa05XVxPQJ3bR7v4bHJhpdGgAkTiKCe+fIuLIZU3dHu8yCz9mMaecIU44AQKlEBPfBqRl1ZTNnbOvKZnRoaqZBFQFAciUiuNf3dGt2PnfGttn5nNb1dDeoIgBIrkQE9/bN/ZrPuWbmFuQefJ7PubZv7m90aQCQOIkI7qGBPu3Yukl9qzp1dHZefas6tWPrJkaVAEAZiRkOODTQR1ADQASJaHEDAKIjuAEgZQhuAEiZSMFtZmvM7A4zGzOzA2b2irgLAwCUF/XNyU9Kutvd/9jMOiQxwBoAGqRicJvZakmbJb1Dktx9TtJcvGUBABYTpaukX9KkpP80sx+b2c1mdl7pQWa2zcxGzWx0cnKy5oUCAAJRgrtd0pWSPuvuV0h6RtKHSw9y913uPujug729kRYqBgCcgyjBfUjSIXe/L/z6DgVBDgBogIp93O7+KzM7aGbPd/dHJL1K0k/iKogFFQBgaVHHcb9X0hfM7CFJl0v6xziKYUEFAKgs0nBAd39A0mDMtZyxoIIkdXe0a2ZuQTtHxml1A0AoUXdOsqACAFSWqOBmQQUAqCxRwc2CCgBQWaKCmwUVAKCyxCykUMCCCgCwtES1uAEAlRHcAJAyBDcApAzBDQApQ3ADQMoQ3ACQMgQ3AKQMwQ0AKUNwA0DKENwAkDIENwCkDMENAClDcANAyhDcAJAyBDcApAzBDQApQ3ADQMoQ3ACQMpGWLjOzn0malpSTtODug3EWNTw2oZ0j4zo4NaP1Pd3avrmf5cwAIFRNi/tqd7+8HqF90+79mpg+oTVdWU1Mn9BNu/dreGwizh8LAKmRuK6SnSPjymZM3R3tMgs+ZzOmnSPjjS4NABIhanC7pG+a2V4z21buADPbZmajZjY6OTl5zgUdnJpRVzZzxraubEaHpmbO+ZwA0EyiBvcr3f1KSddKutHMNpce4O673H3Q3Qd7e3vPuaD1Pd2anc+dsW12Pqd1Pd3nfE4AaCaRgtvdD4efJyR9VdJL4ypo++Z+zedcM3MLcg8+z+dc2zf3x/UjASBVKga3mZ1nZqsKjyX9vqR9cRU0NNCnHVs3qW9Vp47OzqtvVad2bN3EqBIACEUZDvhbkr5qZoXjv+jud8dZ1NBAH0ENAIuoGNzuPi7psjrUAgCIIHHDAQEASyO4ASBlCG4ASBmCGwBShuAGgJQhuAEgZQhuAEgZghsAUobgBoCUIbgBIGUIbgBIGYIbAFKG4AaAlCG4ASBlCG4ASBmCGwBSJsoKOA0zPDahnSPjOjg1o/U93dq+uZ+VcQC0vMS2uIfHJnTT7v2amD6hNV1ZTUyf0E2792t4bKLRpQFAQyU2uHeOjCubMXV3tMss+JzNmHaOjDe6NABoqMQG98GpGXVlM2ds68pmdGhqpkEVAUAyJDa41/d0a3Y+d8a22fmc1vV0N6giAEiGxAb39s39ms+5ZuYW5B58ns+5tm/ub3RpANBQiQ3uoYE+7di6SX2rOnV0dl59qzq1Y+smRpUAaHmRhwOaWUbSqKQn3f318ZV02tBAH0ENACWqaXG/X9KBuAoBAEQTKbjNbJ2k10m6Od5yAACVRG1xf0LS30jKL3aAmW0zs1EzG52cnKxJcQCAs1UMbjN7vaQJd9+71HHuvsvdB919sLe3t2YFAgDOFKXF/UpJW83sZ5Jul3SNmd0aa1UAgEVVDG53/1t3X+fuGyVdL+ked39r7JUBAMpK7DhuAEB5VU3r6u7DkoZjqQQAEAktbgBIGYIbAFKG4AaAlCG4ASBlCG4ASBmCGwBSJtGrvEus9A4ApRLd4maldwA4W6KDm5XeAeBsiQ5uVnoHgLMlOrhZ6R0Azpbo4GaldwA4W6KDm5XeAeBsiR8OyErvAHCmRLe4AQBnI7gBIGUIbgBIGYIbAFKG4AaAlCG4ASBlCG4ASBmCGwBSpmJwm1mnmf3QzB40s/1m9g/1KAwAUF6UOydPSrrG3Y+bWVbS98zsLnffE3NtAIAyKga3u7uk4+GX2fDD4yyqHFbCAYBApD5uM8uY2QOSJiR9y93vK3PMNjMbNbPRycnJmhbJSjgAcFqk4Hb3nLtfLmmdpJea2YvKHLPL3QfdfbC3t7emRbISDgCcVtWoEnd/WtKwpC2xVLMIVsIBgNOijCrpNbM14eMuSa+WNBZ3YcVYCQcATovS4r5I0r1m9pCkHyno4/5avGWdiZVwAOC0KKNKHpJ0RR1qWdTQQJ92KOjrPjQ1o3WMKgHQwhK/Ak4BK+EAQIBb3gEgZQhuAEgZghsAUobgBoCUIbgBIGUIbgBIGYIbAFImNeO4C5jeFUCrS1WLm+ldASBlwc30rgCQsuBmelcASFlwM70rAKQsuJneFQBSFtxDA33asXWT+lZ16ujsvPpWdWrH1k2MKgHQUlI3HJDpXQG0ulS1uAEABDcApA7BDQApk7o+bonb3gG0ttS1uLntHUCrS11wc9s7gFaXuuDmtncAra5icJvZejO718wOmNl+M3t/PQpbDLe9A2h1UVrcC5L+yt1fIOnlkm40sxfGW9biuO0dQKurGNzu/kt3vz98PC3pgKRnx13YYrjtHUCrq2o4oJltlHSFpPvK7NsmaZskbdiwoQalLY7b3gG0sshvTprZSklfkfQBdz9Wut/dd7n7oLsP9vb21rJGAECRSC1uM8sqCO0vuPud8ZYUHTfiAGhFUUaVmKT/kHTA3T8ef0nRcCMOgFYVpavklZLeJukaM3sg/Lgu5roq4kYcAK2qYleJu39PktWhlqocnJrRmq7sGdu4EQdAK0jdnZMF3IgDoFWlNri5EQdAq0ptcHMjDoBWlcr5uAsKIV0YElh4Y5LwBtDMUtvilhgSCKA1pTq4GRIIoBWlOriZmxtAK0p1cDMkEEArSnVwMyQQQCtKdXAXDwn81dFZTU6f1MzcgnaOjPMGJYCmlergloLw3r65X90rsupdtUIXru5kdAmAppb64JYYXQKgtTRFcDO6BEAraYrgLh5dcmx2XuOTxzX2q2kdnZ2nuwRA02mK4C6MLpmcPqHDR2c1l8vLJJ23IkNfN4Cm0xTBXRhdMjOXU95dHZk2PbunS2tXdtLXDaDppHqSqWJDA31a3ZXVhgu6Fay2FqCvG0CzaZrgloK+7onpE1rIuX59/KTmcnll2kwbL+BOSgDNoym6Sgq2b+7X0dl5Pfn0rObDfu6FnOvIM3P0cwNoGk0V3EMDfepduULtbSaXlM20aV1Pl1Z3ZennBtA0mqqrRJKmTy7oeX0rZWY6Njuvyemgy+TQ1KyGxyZYZAFA6jVVi1s6Pab72Oy8Dh+d1ULOZQqWqWdoIIBm0HTBXRjT/dT0CZmknOc1n3PN510T0yf0sbsONLpEAFiWisFtZreY2YSZ7atHQctVGNPtHrwxmctLmTZTts2Uz7semzxOqxtAqkVpcX9O0paY66ipoYE+XbmhRx3ZNnVk2pTNtCmXD1rdubz0vtt/THgDSK2Kwe3uI5J+U4daaqrQZeJyLeRyms+73KVsm/TM3AL93QBSq2Z93Ga2zcxGzWx0cnKyVqc9Z0MDfbq0b6XazLSQD96c7Mi0ySW5S4efnqXlDSCVahbc7r7L3QfdfbC3t7dWp12WD20ZUN/qTmXaTB3tprzntZCX5FLOXcdOLGj7rXv1qW8/2uhSASCyphtVUqzwRmV3R0YLeVfeg184r6DVLUlzC3l9/NuP6dpPjND6BpAKTR3cUhDen7r+Cj17TbfazOR25v4wv/XIU9O64fOjBDiAxIsyHPA2ST+Q9HwzO2RmN8RfVm0Vt7zddeqGnGJ5l/J514FfTdN9AiDRzAt9BjU0ODjoo6OjNT/vcg2PTWj7rXuVd9d8bvHf2xS0xF9w4Sp9aMsAt8kDiJ2Z7XX3wSjHNn1XSbGhgT7dOHSx2sL5uktb3QWFSKf1DSCJWqrFXTA8NqGP3XVAj00el3vQTbKYQuu7uyOj7o6MLulbpe2b+2mFA6ipalrcTTc7YBRDA30aGug7FeBjTx1f9NhCps/M5TQ7l9PUM0d0/y+mdOPQxXrfqy+tT8EAUKQlW9ylPvXtR/Xp4Z/q5EK+qu/Ltpme17eSfnAAy1ZNi5vgDkVpfZdDVwqAWiC4l6G49V148zLqFbKiB6tWtGvTb59PiAOIhD7uZXjfqy/Vi9etOfXmZZuk+Yg9KF704NiJBe0ZP6I9TxxRZ3swQ2FHexstcgDLRot7CcNjE9o5Mq59T07p+Ml85JZ3OcVDD13BHOHtbab+tefRRw6ArpI4FPrAH3nquFyn+7aXoy08h0vqyAQLHLcZYQ60IoI7RoUAf3zyuBbyUsakJW7CrFq5MHd3dbTz5ifQzAjuOintSpGW3wqvxML/dLYHN70u5J1gB5oAwd0A5UK8IO4wL0awA+lEcDdYIcQfe+qY5nKuuYWcTi4E17meIV7OYsHeZnTLAI1EcCdQcZg/M5fTyYW88h70kee98YFeyiS1hf33bZI6OzLKWBD0hD1QewR3SqQtzCspDfv2kjdXCX5gcQR3ipWG+ULelcv5WSGetlCPotrg58UAzYTgbkKl/ebufiqwWiXYz8Vy/wqoZj9j8LEcBHcLItgbb7Ex+Mt5MajFfv7ySAeCG2dZKthL/+dvM9NceFdRm6TqJrtFEtXzL4+oLya8QJ2J4MayFe4QfeLIjHL5fMX/kQh7JEXUIa9xvJgs54WD4EZDlIZ9NU96gh/NIGNS3+oVymYy2rF1U1XhzbSuaIjCknDnqtrg58UASWIKhvEem13Qhee3a+fIeGzdNZGC28y2SPqkpIykm939Y7FUg5a23OAvZzl/BVSzP+1j8FEDJrlLc7m8urIZHZqaie1HVQxuM8tI+rSk10g6JOlHZrbb3X8SW1VAjcTxYrCYcmPwkzCqhL886sSDVndHpk2z8zmt6+mO7UdFaXG/VNLj7j4uSWZ2u6Q3SCK4gSL1fJGoVr3+8mjUqJIkDHl1BX9xre5q13zOtX1zf2w/K0pwP1vSwaKvD0l6WelBZrZN0jZJ2rBhQ02KA1AbSX5RqZVqhrzGOapk47NWxj4cMUpwW5ltZ72YufsuSbukYFTJMusCgKq0wotTQVuEYw5JWl/09TpJh+MpBwBQSZTg/pGkS8zsuWbWIel6SbvjLQsAsJiKXSXuvmBm75H0DQXDAW9x9/2xVwYAKCvSOG53/7qkr8dcCwAggihdJQCABIllrhIzm5T085qfuDbWSvp1o4tYQtLrk5JfI/UtD/Ut37nU+Bx3741yYCzBnWRmNhp1IpdGSHp9UvJrpL7lob7li7tGukoAIGUIbgBImVYM7l2NLqCCpNcnJb9G6lse6lu+WGtsuT5uAEi7VmxxA0CqEdwAkDJNFdxmtsXMHjGzx83sw2X2/6WZ/cTMHjKz75jZc4r25czsgfAjlrlYItT3DjObLKrjXUX73m5mj4Ufb29Qff9aVNujZvZ00b56XL9bzGzCzPYtst/M7FNh/Q+Z2ZVF++px/SrV95awrofM7PtmdlnRvp+Z2cPh9YtlwdYI9Q2Z2dGif8ebivYt+dyoU30fLKptX/icuyDcV4/rt97M7jWzA2a238zeX+aY+jwH3b0pPhTMo/JTSf2SOiQ9KOmFJcdcLak7fPznkr5UtO94Aup7h6R/K/O9F0gaDz/3hI976l1fyfHvVTBvTV2uX/gzNku6UtK+RfZfJ+kuBVMRv1zSffW6fhHru6rwcyVdW6gv/PpnktY2+PoNSfracp8bcdVXcuwfSLqnztfvIklXho9XSXq0zP/DdXkONlOL+9RKPe4+J6mwUs8p7n6vuxcWgtujYIraxNS3hNdK+pa7/8bdpyR9S9KWBtf3p5Juq3ENS3L3EUm/WeKQN0j6vAf2SFpjZhepPtevYn3u/v3w50v1f/5FuX6LWc5zN7Iq62vE8++X7n5/+Hha0gEFC80Uq8tzsJmCu9xKPaUXtdgNCl4ZCzrNbNTM9pjZHzawvj8K/8S6w8wK86BX+7vFWZ/CLqbnSrqnaHPc1y+KxX6Hely/apU+/1zSN81srwWrSTXKK8zsQTO7y8w2hdsSdf3MrFtB6H2laHNdr5+ZbZR0haT7SnbV5TkYaXbAlIi0Uo8kmdlbJQ1K+r2izRvc/bCZ9Uu6x8wedvef1rm+/5F0m7ufNLN3S/ovSddE/N561FdwvaQ73D1XtC3u6xfFYr9DPa5fZGZ2tYLg/t2iza8Mr1+fpG+Z2VjYAq2n+xXMl3HczK6T9N+SLlHCrp+CbpL/c/fi1nndrp+ZrVTwovEBdz9WurvMt9T8OdhMLe5IK/WY2aslfUTSVnc/Wdju7ofDz+OShhW8mta1Pnc/UlTTv0t6SdTvrUd9Ra5XyZ+pdbh+USz2OyRmFScze7GkmyW9wd2PFLYXXb8JSV9V0D1RV+5+zN2Ph4+/LilrZmuVoOsXWur5F+v1M7OsgtD+grvfWeaQ+jwH4+zMr+eHgr8exhX8CV94A2VTyTFXKHiT5ZKS7T2SVoSP10p6TDV+8yVifRcVPX6jpD1++o2NJ8I6e8LHF9S7vvC45yt4I8jqef2KftZGLf7m2ut05htDP6zX9YtY3wZJj0u6qmT7eZJWFT3+vqQtDajvwsK/q4Lg+0V4LSM9N+KuL9x/voJ+8PPqff3Ca/F5SZ9Y4pi6PAdrfuEb+aHgHd1HFYTzR8JtOxS0riXp25KekvRA+LE73H6VpIfDJ+TDkm5oUH3/JGl/WMe9kgaKvved4f/0j0v6s0bUF37995I+VvJ99bp+t0n6paR5BS2YGyS9W9K7w/0m6dNh/Q9LGqzz9atU382Spoqef6Ph9v7w2j0Y/vt/pEH1vafo+bdHRS8w5Z4b9a4vPOYdkm4v+b56Xb/fVdC98VDRv+F1jXgOcss7AKRMM/VxA0BLILgBIGUIbgBIGYIbAFKG4AaAJVSa/OoczvfP4SRZ+8zsT87lHAQ3ACztc6rR3DZm9joFE2ldLullkj5oZqurPQ/BDQBL8DKTX5nZxWZ2dzg3ynfNbCDi6V4o6X/dfcHdn1Ew9rzqFwWCGwCqt0vSe939JZL+WtJnIn7fg5KuNbPucDqBq3XmrfCRNNMkUwAQu3CSqaskfdns1NxRK8J9b1Jwt3GpJ939te7+TTP7HQW35U9K+oGkhapr4M5JAFhaOI3r19z9RWGf9CPuflENzvtFSbd6MKlXZHSVAEAVPJjK9Qkze7N0armyyyp8m8JjM2b2rPDxiyW9WNI3q62BFjcALMHMblOwrNtaBZPUfVTBIiKfVbCcWVbBxFflukhKz9WpYN5zSTqmYHKqB6quieAGgHShqwQAUobgBoCUIbgBIGUIbgBIGYIbAFKG4AaAlCG4ASBl/h88u2J+65EXFQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "R=[i*Ao/10 for i in range(10,200)]\n",
    "Y=[Coul_Buck(r, (r,0,0), Cl,Cl,NaCl)[0] for r in R]\n",
    "sns.regplot(x=R,y=Y,fit_reg=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x120236e80>"
      ]
     },
     "execution_count": 50,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEQCAYAAABLMTQcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGqBJREFUeJzt3X1wXOd13/HvDyBIACIpQiZpqgRZmrESxvL4hWVkNdawjNwkspKRklTpiI2dKlVCKmnivqZKOhk5ZfuPJjNpxq1lk000SuJacuI4LquRpaaRNVRbUzEp2Y5k0TGHUiRQlgBKfAEFgAAWp3/sAlwsF9hd4O7uvXd/nxkOsXsv7z66FA4PzvM85yoiMDOzfOlq9wDMzCx5Du5mZjnk4G5mlkMO7mZmOeTgbmaWQw7uZmY51NbgLulBScOSnk/oeo9LOifp0Yr3PyLpWUnfkPR/JL07ic8zM0urdmfuDwG3JHi93wE+XuX9zwA/FxEfAD4P/FaCn2lmljptDe4RcQR4q/w9Sd9XysCPS3pa0o4GrveXwGi1Q8Da0tdXA68tdcxmZlmwot0DqOIQcE9EfFfSh4AHgJuXec1fBB6TNA5cAG5c5vXMzFItVcFd0mrgh4E/lTT79qrSsZ8BDlT5Y6cj4sdrXPpfAbdGxDOSfh34XYoB38wsl1IV3CmWic6VauPzRMSXgC81ekFJG4D3R8Qzpbe+ADy+rFGamaVcuydU54mIC8BLkn4WQEXvX+ZlzwJXS/r+0usfBV5c5jXNzFJN7ewKKelhYA+wHngD+CTwJMXVLdcCPcAjEVGtHFPtek8DO4DVwJvA3RHxhKSfpljSmaEY7P9ZRJxK9r/GzCw92hrczcysOVJVljEzs2S0bUJ1/fr1sW3btnZ9vJlZJh0/fvxMRGyodV7bgvu2bds4duxYuz7ezCyTJP1tPee5LGNmlkMO7mZmOeTgbmaWQw7uZmY55OBuZpZDaest01RPnRjm4JFTvHp2jC0D/ezfvZ09Oza2e1hmZonrmMz9qRPD3Hf4BYZHJ1jX18Pw6AT3HX6Bp04Mt3toZmaJ65jgfvDIKXq6Rf/KFUjF33u6xcEjbjFjZvmT67JMeRlmZPQSm9aumnd8ujDDs6+c5ab7n3SZxsxyJbeZe2UZRoLT5yYYnZgC4ML4FKfPTSBwmcbMcie3wb2yDPPONb0AvH5+gojgjdEJADZd3esyjZnlTm6D+6tnx+jr6Z57vbavh83regng/PgUEbB5XS9renvmzunr6Wbo7FgbRmtmlqzcBvctA/2MTxXmvbeiu4udWwd4+t6b2bl1gBXdl//zL4xPcXLkIsOjl9h76KjLM2aWabkN7vt3b2eqEIxNThNR/H2qEOzfvf2K4xfGJzl9bpzpQrBp7SrX380s83Ib3Pfs2MiB265n45pezo9PsXFNLwduu35uNUz58dcvXGJFlxgc6GNt30rX380s83K9FHLPjo2LLm2cPX7T/U+WVtRo7pjr72aWZbnN3BtRrT4/PlVgcKC/TSMyM1uemsFd0oOShiU9X+O8H5JUkHRHcsNrjVr1eTOzrKknc38IuGWxEyR1A/cDTyQwpparVp+/Y+dmDh45xU33P+nVM2aWOTVr7hFxRNK2Gqf9GvBnwA8lMKa2KK/Pz+5u7enWvN2rB0rnmZml3bInVCVtBn4auJkawV3SPmAfwNatW5f70VdIqqVv+e5WgP6VKxibnObgkVMO7maWCUlMqP4ecG9EFGqdGBGHImJXROzasGFDAh99WZItfSt3t4JXz5hZtiQR3HcBj0h6GbgDeEDSTyVw3YYk2dLXq2fMLOuWHdwj4l0RsS0itgFfBH4lIr687JE1KMlsu9rqmfPjU5wbm/QEq5llQj1LIR8Gvgb8gKQhSXdLukfSPc0fXv2SzLYrV8/0dAkBk4UZtwc2s0yoZ7XM3novFhF3LWs0y7B/93buO/wCY5PT9PV0Mz5VWNZa9fLVM3sPHWVqJjzBamaZkZsdqrV6ySyHJ1jNLGty1VumVi+Zpdoy0M/w6MRc5g6eYDWzdMtN5t5Mbk9gZlnj4F4Htycws6zJVVmmmdyewMyyxJn7EiS5YcrMrBkc3JfAq2fMLO0c3JfA7QnMLO0c3JfA7QnMLO0c3JfA7QnMLO28WmaJ3J7AzNLMmXsCPMFqZmnj4J4AT7CaWdo4uCegcoJ1ZHSCobPjfHd41JOrZtYWDu4JKJ9gff38OGfHprjmqh42re315KqZtYWDe0L27NjIw/tu5Lp3rmVwoI/1q3u9e9XM2sbBPWGeXDWzNHBwT5gnV80sDRzcE+bdq2aWBg7uCfPuVTNLA+9QbQLvXjWzdnPm3mSeYDWzdnBwbzJPsJpZO2Q6uD91Ypi9h46meqLSE6xm1g6ZDe6zzzEdHp1I9USlJ1jNrB0yG9yz9BzT2d2rT997MwNXrWJtX08mxm1m2ZXZ4J7VicqsjtvMsiWzwT2rE5XVxn3m4iXOj0+5Bm9miclscK82UTlVCPbv3t7uoS2qWnvgkYuTXLWq2zV4M0tMZoN75UTlxjW9HLjt+tRvDKoc99hkgY1rVrqLpJklquYOVUkPAj8JDEfEe6sc/zng3tLLi8AvR8Q3Ex3lAsp3gmZJ+bhvuv9J1vX1zDvuGryZLVc9mftDwC2LHH8J+AcR8T7gPwKHEhhXx8jq3IGZpVvN4B4RR4C3Fjn+/yLibOnlUWAwobF1BG9yMrNmSLrmfjfwlYUOSton6ZikYyMjIwl/dDZ5k5OZNYMiovZJ0jbg0Wo197JzfgR4ALgpIt6sdc1du3bFsWPH6h9ph9h76CjDoxNzXSQBxian2biml4f33djGkZlZGkg6HhG7ap2XSOYu6X3A7wO31xPYbWHe5GRmSVh2cJe0FfgS8PGI+JvlD6mzVU6wXhif4uTIRYZHL7n+bmZ1qxncJT0MfA34AUlDku6WdI+ke0qn3Ae8A3hA0jckudayDOUTrBfGJzl9bpzpQrBp7SrX382sbnXV3JvBNfeFPXVimINHTvHsK2cRsOnqXtb0FtfCu/5u1tlaWnO3ZM12kdywZhXv3rh6LrCD6+9mVh8H9xRzkzEzWyoH9xRzkzEzWyoH9xRzkzEzW6qajcOsvdxkzMyWwpl7hrjJmJnVy8E9Q9xkzMzq5eCeIW4yZmb1cs09Y8pr8HsPHWVqJuaajPWvXMHY5DQHj5zK5ENMzCw5ztwzzE3GzGwhDu4Z5k1OZrYQB/cM8yYnM1uIg3uGeZOTmS3EE6oZ501OZlaNM/cccQ3ezGY5uOeIa/BmNsvBPUdcgzezWa6554xr8GYGztxzzQ/bNutcDu455odtm3UuB/ccK6/Bv37hEiu6xOBAH2v7Vrr+bpZzDu4554dtm3UmB/cO4TXwZp3Fwb1DeA28WWdxcO8QXgNv1lm8zr2DeA28Wedw5t6hXIM3yzcH9w7lGrxZvjm4dyjX4M3yzTX3DuYavFl+1czcJT0oaVjS8wscl6RPSTop6VuSdiY/TGs21+DN8qWessxDwC2LHP8ocF3p1z7gM8sflrWaa/Bm+VIzuEfEEeCtRU65HfijKDoKrJN0bVIDtNZwDd4sX5KouW8GXi17PVR673uVJ0raRzG7Z+vWrQl8tCWpVg1+ujDDs6+c5ab7n2TLQD/7d2+fO9/M0iWJ1TKq8l5UOzEiDkXErojYtWHDhgQ+2pqlWi/40+cmELhMY5YBSQT3IWBL2etB4LUErmttVFmDf2N0AoBNV7tMY5YFSQT3w8DPl1bN3Aicj4grSjKWLZU1+AjYvK7XLYPNMqJmzV3Sw8AeYL2kIeCTQA9ARHwWeAy4FTgJjAG/0KzBWmuV1+D3HjrKcCl7h2KZ5o3RCSKKx1x/N0sXRVQtjzfdrl274tixY235bGvcUyeGue/wC/R0i+nCDKfPFQP95nW9rOjuYqoQHLjtegd4syaTdDwidtU6z+0HrC5+ZJ9Ztjhzt4bNLpOULi+UujA+yesXLrFhzSovkzRrImfu1jReJmmWfg7u1jAvkzRLPwd3a5iXSZqln2vutmyzyyT7V15eWTsyOsHYZIG1fT2uwZslyDV3axl3lDRLHwd3WzZ3lDRLH5dlLHFeKmnWPPWWZTL3mL2nTgxz8MgpXj075iCRUlsG+ufV4GeXSq7o0rwyzQHw351Zk2SqLDO7BX54dMK13BTzUkmz9stUcD945BQ93cXg4CCRXvUslSx/8Iefz2qWvEyVZV49O3bF04G8njqdanWUdJnGrLkylblXbnsHGJ8qMDjQ36YRWT1cpjFrvUwF98ogMTY5zVQh2L97e7uHZovwjlaz1svcUsjZ1TJDZ8cY9GqZTKrc0Vr+4I+dWwf8d2q2iHqXQmYuuFv2+cEfZkvn9gOWWn7wh1nzOXO3tvJuVrPGOHO3TPCDP8yaw8Hd2srLJM2aw8Hd2sq7Wc2aI1M7VC2fvJvVLHnO3C1VXKYxS4aDu6WKyzRmyXBZxlLHZRqz5XPmbqlWT5lmcrrAJx55zpm8WRkHd0u1WmWaC+NTvPn2JG9PTntdvFkZl2Us9RYr05y5eAmA3hXdc5n82OQ0B4+ccpnGOpozd8uUyjLNpekZCNiwZhVQzORfPz/BX738lks01tEc3C1TKss0/Su7Wb9mJWt6e7gwPsVr58eZLMywqlsu0VhHqyu4S7pF0ncknZT0G1WOb5X0VUnPSfqWpFuTH6pZ0Z4dG3l43408fe/NfOrOD9LT3c3Y5PRciUaIjWt7vSbeOlrN4C6pG/g08FHgPcBeSe+pOO23gD+JiA8CdwIPJD1Qs2rKM/mJ6Rl6usTfKZtw9Zp461T1ZO43ACcj4lRETAKPALdXnBPA2tLXVwOvJTdEs8XNZvI3bLuGa9f1zVtJ4w6T1qnqCe6bgVfLXg+V3iv328DHJA0BjwG/Vu1CkvZJOibp2MjIyBKGa7Ywr4k3u6ye4K4q71U+4WMv8FBEDAK3An8s6YprR8ShiNgVEbs2bNjQ+GjNFuE18WaX1bPOfQjYUvZ6kCvLLncDtwBExNck9QLrAX/XWEt5TbxZUT2Z+9eB6yS9S9JKihOmhyvOeQX4CICkHwR6AdddrK1qrYkHT7haftUM7hExDfwq8ATwIsVVMS9IOiDpttJp/wb4JUnfBB4G7op2PZzVrGSxNfHgCVfLNz8g2zrGUyeGue/wC/R0i76ebk6OXGS6EAwOXF5hMzI6wdhkgbV9PX44t6WSH5BtVsETrtZJ3DjMOoonXK1TOHO3juUmZJZnDu7WsdyEzPLME6pmJeUTrq+fn2CyMIO43KvGk62WBp5QNWvQYk3IPNlqWeMJVbMysxOus5Ot/SuL3yLVJltHRif4xCPPOZO3VHLmblZFPZOtzuQtzZy5m1WxZ8dGDgAHj5xi6OwY/Su7uWpV99yaeGfylnYO7mYLKF8TPzvZOjY5TV9PN5emZxBXZvIzEWy9pn8ukz9Quo5Zq7ksY1aHWn1qqmXyfsSftZMzd7M6NZLJw/yOky7TWKs5czdbAnectLRz5m62RItl8pWP+JsuBMOjE+z/3HF2bh1wFm9N58zdLAGLdZycbWUwMxMUZmacxVtLOHM3S8hCHSfPXLxEFwLBqu4uL5u0lnDmbtYE5ZugJgszBEGUNkF5A5S1gjN3syaYvwlqHFGsv6/p7eHUyEXAG6CsuRzczZpktkwzO9na3aW5Vga1NkD92y9+kw2rVzF6adrB3pbEZRmzJmt0A9R0ITg3NsVLZ9522caWzJm7WQs0sgHqzMVLdAkKES7b2JI5czdrsVqZ/GRhBoCV3cVvT0/A2lI4czdrg8Uy+e6uYmmmPJOH+WUbb4iyWpy5m7VZZSa/7Zp+Bvp75k3AUraM0huirB7O3M1SoDyTh2I2X62X/KmRi94QZXXxA7LNUq78wd2vvDWGAEoP7o6gmMlH8IOb1jI+VeD8+JSXUeaYH5BtlhPlZZsuiS5dfnC3l1HaQlyWMcuAejdEeRmlzXLmbpYhXkZp9XLmbpYxy1lG6Uy+c9SVuUu6RdJ3JJ2U9BsLnPOPJX1b0guSPp/sMM2smkaWUcKVmfxLZy6y/3PH2fWf/oK9h446o8+Rmpm7pG7g08CPAkPA1yUdjohvl51zHfCbwIcj4qwkpwFmLVLvMkqYn8mPTkzz5tuTAIxdmp4r2RwoXdOyrZ6yzA3AyYg4BSDpEeB24Ntl5/wS8OmIOAsQEf7n36xN6u1jMzJafIiIumBqJlyyyZl6gvtm4NWy10PAhyrO+X4ASf8X6AZ+OyIer7yQpH3APoCtW7cuZbxm1oD5feXnZ/Knz43TLRFRnIB16+F8qSe4q8p7lTufVgDXAXuAQeBpSe+NiHPz/lDEIeAQFDcxNTxaM2vYQpn8yu4uJgszCM1l8nDlmvmLE9O8e+Nql20ypp7gPgRsKXs9CLxW5ZyjETEFvCTpOxSD/dcTGaWZJaI8kz8/Nsn0THDNVT2sXrVi7olRXjOfD/UE968D10l6F3AauBP4JxXnfBnYCzwkaT3FMs2pJAdqZsmozOQXmnwtZvVXrpl32SYbagb3iJiW9KvAExTr6Q9GxAuSDgDHIuJw6diPSfo2UAB+PSLebObAzWz5ltt62GWb9HLjMDObU57JX7WymzffnmRtXw99Pd2ceH0UAZsH+uY6VE4VZghgx6a1XBif4o3RCSJwn/kmqrdxmHeomtmcRtbMl5dtZvvMC5iJcMkmBRzczWxB9ZZtZtfMz/aZr1aycbBvLTcOM7O6LNbqYLIwQxBEqdVB5UobtyJuPWfuZla3hcs2xZLMpqt75zZIla+08bLK1nNwN7MlW6jPfOVKGy+rbD0HdzNbtso2B9uu6efNtycXDPa1llXOdqtc07uC6zaucaBfAgd3M0vEYittKoP9Yk+Qqtat0ll94xzczawplrqs8szF+d0qvfJmaRzczawl6l1WWdmtcqGVN94ZuzgHdzNrucVq9JXdKr3yZmkc3M2sLRYq21R2q/TKm6VxcDezVFioW6VX3iyNg7uZpY5X3iyfg7uZpZ5X3jTOwd3MMqeZK2/yEuwd3M0s05JceZOnYO/gbmaZl9TKm1rBPkuTsw7uZpY7S115s1iwz9rkrIO7meVaIytvFgv2WZucdXA3s46y1GCftclZB3cz62j1BvusTc46uJuZlWnV5GyzG545uJuZLaIZk7OzDc/GJqc5eOSUg7uZWTslNTk7q6+nm6GzY00Zq4O7mdkSLSfYA4xPFRgc6G/K2BzczcwS0kiwH58qMFUI9u/e3pSxOLibmTXJYsF+0KtlzMzyoTLYN1NX7VPMzCxrHNzNzHLIwd3MLIcc3M3McsjB3cwshxQR7flgaQT424Qutx44k9C1kuRxNcbjaozH1Zi8jOvvRsSGWie1LbgnSdKxiNjV7nFU8rga43E1xuNqTKeNy2UZM7MccnA3M8uhvAT3Q+0ewAI8rsZ4XI3xuBrTUePKRc3dzMzmy0vmbmZmZRzczcxyKFPBXdKDkoYlPb/AcUn6lKSTkr4laWdKxrVH0nlJ3yj9uq8FY9oi6auSXpT0gqR/UeWclt+vOsfVjvvVK+mvJH2zNK7/UOWcVZK+ULpfz0jalpJx3SVppOx+/WKzx1X22d2SnpP0aJVjLb9fdY6rLfdL0suS/rr0mceqHE/2+zEiMvML2A3sBJ5f4PitwFcAATcCz6RkXHuAR1t8r64Fdpa+XgP8DfCedt+vOsfVjvslYHXp6x7gGeDGinN+Bfhs6es7gS+kZFx3Af+1lfer7LP/NfD5an9f7bhfdY6rLfcLeBlYv8jxRL8fM5W5R8QR4K1FTrkd+KMoOgqsk3RtCsbVchHxvYh4tvT1KPAisLnitJbfrzrH1XKle3Cx9LKn9KtytcHtwB+Wvv4i8BFJSsG42kLSIPATwO8vcErL71ed40qrRL8fMxXc67AZeLXs9RApCBwlf7/0o/VXJF3fyg8u/Tj8QYpZX7m23q9FxgVtuF+lH+W/AQwDfxERC96viJgGzgPvSMG4AP5R6Uf5L0ra0uwxlfwe8O+AmQWOt+V+1TEuaM/9CuB/STouaV+V44l+P+YtuFfLCtKQ5TxLsR/E+4H/Any5VR8saTXwZ8C/jIgLlYer/JGW3K8a42rL/YqIQkR8ABgEbpD03opT2nK/6hjX/wS2RcT7gP/N5Wy5aST9JDAcEccXO63Ke029X3WOq+X3q+TDEbET+CjwzyXtrjie6P3KW3AfAsr/FR4EXmvTWOZExIXZH60j4jGgR9L6Zn+upB6KAfS/R8SXqpzSlvtVa1ztul9ln38OeAq4peLQ3P2StAK4mhaW4xYaV0S8GRGXSi//G/D3WjCcDwO3SXoZeAS4WdLnKs5px/2qOa423S8i4rXS78PAnwM3VJyS6Pdj3oL7YeDnS7PONwLnI+J77R6UpE2ztUZJN1C87282+TMF/AHwYkT87gKntfx+1TOuNt2vDZLWlb7uA/4hcKLitMPAPy19fQfwZJRmwto5roq67G0U5zGaKiJ+MyIGI2IbxcnSJyPiYxWntfx+1TOudtwvSVdJWjP7NfBjQOXqukS/HzP1gGxJD1NcSbFe0hDwSYoTTETEZ4HHKM44nwTGgF9IybjuAH5Z0jQwDtzZ7P/JKWYwHwf+ulSvBfj3wNaycbXjftUzrnbcr2uBP5TUTfEfkz+JiEclHQCORcRhiv8o/bGkkxQz0DubPKZ6x/UJSbcB06Vx3dWCcVWVgvtVz7jacb/eCfx5KWdZAXw+Ih6XdA805/vR7QfMzHIob2UZMzPDwd3MLJcc3M3McsjB3cwshxzczcwaoBqNApdwvcclnVNFkzNJ71Kx4dp3VWzAtrKR6zq4m5k15iGu3OC2HL9DcXlwpfuB/xwR1wFngbsbuaiDu5lZA6o1CpT0faUM/LikpyXtaOB6fwmMVlxPwM0UG65BsUXCTzUyzkxtYjIzS6lDwD0R8V1JHwIeoBicl+odwLlSwzVYQhMxB3czs2UoNcH7YeBPyzoaryod+xngQJU/djoifnyxy1Z5r6Edpw7uZmbL00Uxy/5A5YFSY7xqTftqOUOxn/uKUvbecBMx19zNzJah1LL6JUk/C3OPy3v/Mq8ZwFcp9lmCYgO2/9HINdxbxsysAeWNAoE3KDYKfBL4DMVGbz3AIxFRrRxT7XpPAzuA1RS7n94dEU9I2k6xbfE1wHPAx8paFde+roO7mVn+uCxjZpZDDu5mZjnk4G5mlkMO7mZmOeTgbmaWQw7uZmY55OBuZpZD/x/J/hjX1wTEBgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "R=[i*Ao/20 for i in range(19,100)]\n",
    "Y=[Coul_Buck(r, (r,0,0), Na,Na,NaCl)[0] for r in R]\n",
    "sns.regplot(x=R,y=Y,fit_reg=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x11fcbdda0>"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX8AAAEQCAYAAABY5nToAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAHmJJREFUeJzt3X90XGd95/H3xyPJkiIltmMp8cZ2HUM4PnEPv6qTBmiNF8I2wB6bsoQmW9qwm9TOUg7dZWFJF9Y96549G+BsYXcbWrvAaWCXJiELRMsG0vzAa7qHpJEDpZiYxhUpduJGSqLYciVFo/F3/5ixMhpfWT/mjmY09/M6x2dm7jy69+sb5zOPnnnucxURmJlZtqyodwFmZrb0HP5mZhnk8DczyyCHv5lZBjn8zcwyyOFvZpZBDR/+kr4oaUjSj1La37clvSjpmxXb3yrpcUk/kPQXkl6ZxvHMzBpRw4c/8KfAtSnu79PAbyRs/yPg1yPitcBXgE+keEwzs4bS8OEfEQeBF8q3SXpFqQd/SNJ3JW1ZwP4eAkaT3gIuLD2/CHhmsTWbmTW6lnoXsEj7gVsi4klJvwh8DnhLlfu8GbhP0jhwCri6yv2ZmTWsZRf+krqANwJflXR288rSe+8G9ib82NMR8Stz7PrfAO+IiEclfRT4A4ofCGZmTWfZhT/FoaoXS2PzM0TE14CvLXSHknqA10TEo6VNdwHfrqpKM7MG1vBj/pUi4hTwU0nXAajoNVXudgS4SNKrSq/fBjxR5T7NzBqWGn1VT0l/BmwH1gLPAr8HPExxds46oBW4MyKShnuS9vddYAvQBTwP3BQR90v6VYpDRmcofhj8y4gYTPdvY2bWGBo+/M3MLH3LbtjHzMyq17Bf+K5duzY2bdpU7zLMzJaVQ4cOPRcRPXO1a9jw37RpEwMDA/Uuw8xsWZH0d/Np52EfM7MMcvibmWWQw9/MLIMc/mZmGeTwNzPLoIad7bNYB44Mse/gIMdGxtiwupPd2zazfUtvvcsyM2soTdXzP3BkiD39hxkanWBVRytDoxPs6T/MgSND9S7NzKyhpBL+kq6V9BNJRyXdOkub90r6saTDkr6SxnEr7Ts4SGtOdLa1IBUfW3Ni30Ev0WNmVq7qYR9JOeB2iithHgcek9QfET8ua3MF8LvAmyJiRFJNxmGOjYyxqqN1xraO1hzHR8ZqcTgzs2UrjZ7/VcDRiBiMiEngTmBnRZvfAm6PiBGAiKjJOMyG1Z2M5wszto3nC6xf3VmLw5mZLVtphP9lwLGy18dL28q9CniVpP8n6RFJiTdkl7RL0oCkgeHh4QUXsnvbZvKFYGxyiojiY74Q7N62ecH7MjNrZmmEvxK2Va4T3QJcQXFd/huAz0tadc4PReyPiL6I6OvpmXNdonNs39LL3h1b6e1u5+R4nt7udvbu2OrZPmZmFdKY6nkc2FD2ej3wTEKbRyIiT/EuXD+h+GHwWArHn2H7ll6HvZnZHNLo+T8GXCHpckltwPVAf0WbbwD/GEDSWorDQJ6CY2ZWJ1WHf0RMAR8E7qd439u7I+KwpL2SdpSa3Q88L+nHwHeAj0bE89Ue28zMFqdhb+PY19cXXs/fzGxhJB2KiL652jXVFb5mZjY/Dn8zswxy+JuZZZDD38wsgxz+ZmYZ5PA3M8sgh7+ZWQY5/M3MMsjhb2aWQQ5/M7MMcvibmWWQw9/MLIMc/mZmGeTwNzPLIIe/mVkGOfzNzDLI4W9mlkEOfzOzDEol/CVdK+knko5KuvU87d4jKSTNeYsxMzOrnarDX1IOuB14O3AlcIOkKxPadQMfAh6t9phmZladNHr+VwFHI2IwIiaBO4GdCe1+H/gUMJHCMc3MrApphP9lwLGy18dL26ZJeh2wISK+eb4dSdolaUDSwPDwcAqlmZlZkjTCXwnbYvpNaQXwGeDfzrWjiNgfEX0R0dfT05NCaWZmliSN8D8ObCh7vR54pux1N/DzwAFJTwFXA/3+0tfMrH7SCP/HgCskXS6pDbge6D/7ZkScjIi1EbEpIjYBjwA7ImIghWObmdkiVB3+ETEFfBC4H3gCuDsiDkvaK2lHtfs3M7P0taSxk4i4D7ivYtueWdpuT+OYZma2eL7C18wsg1Lp+TeyA0eG2HdwkGMjY2xY3cnubZvZvqW33mWZmdVVU/f8DxwZYk//YYZGJ1jV0crQ6AR7+g9z4MhQvUszM6urpg7/fQcHac2JzrYWpOJja07sOzhY79LMzOqqqcP/2MgYHa25Gds6WnMcHxmrU0VmZo2hqcN/w+pOxvOFGdvG8wXWr+6sU0VmZo2hqcN/97bN5AvB2OQUEcXHfCHYvW1zvUszM6urpg7/7Vt62btjK73d7Zwcz9Pb3c7eHVs928fMMq/pp3pu39LrsDczq9DUPX8zM0vm8DczyyCHv5lZBjn8zcwyyOFvZpZBDn8zswxy+JuZZZDD38wsgxz+ZmYZlEr4S7pW0k8kHZV0a8L7H5b0Y0k/lPSQpJ9L47hmZrY4VYe/pBxwO/B24ErgBklXVjT7PtAXEa8G7gE+Ve1xzcxs8dLo+V8FHI2IwYiYBO4EdpY3iIjvRMTZRfQfAdancFwzM1ukNML/MuBY2evjpW2zuQn4VtIbknZJGpA0MDw8nEJpZmaWJI3wV8K2SGwovQ/oAz6d9H5E7I+Ivojo6+npSaE0MzNLksaSzseBDWWv1wPPVDaSdA3wceDNEfFSCsc1M7NFSqPn/xhwhaTLJbUB1wP95Q0kvQ7YB+yIiKEUjmlmZlWoOvwjYgr4IHA/8ARwd0QclrRX0o5Ss08DXcBXJf1AUv8suzMzsyWQyp28IuI+4L6KbXvKnl+TxnHMzCwdvsLXzCyDHP5mZhnk8DczyyCHv5lZBjn8zcwyyOFvZpZBqUz1XG4OHBli38FBjo2MsWF1J7u3bWb7lt56l2VmtmQy1/M/cGSIPf2HGRqdYFVHK0OjE+zpP8yBI77w2MyyI3Phv+/gIK050dnWglR8bM2JfQcH612amdmSyVz4HxsZo6M1N2NbR2uO4yNjs/yEmVnzyVz4b1jdyXi+MGPbeL7A+tWddarIzGzpZS78d2/bTL4QjE1OEVF8zBeC3ds217s0M7Mlk7nw376ll707ttLb3c7J8Ty93e3s3bHVs33MLFMyOdVz+5Zeh72ZZVrmev5mZubwNzPLJIe/mVkGOfzNzDIolfCXdK2kn0g6KunWhPdXSrqr9P6jkjalcVwzM1ucqsNfUg64HXg7cCVwg6QrK5rdBIxExCuBzwCfrPa4Zma2eGn0/K8CjkbEYERMAncCOyva7ATuKD2/B3irJKVwbDMzW4Q0wv8y4FjZ6+OlbYltImIKOAlcXLkjSbskDUgaGB4eTqE0MzNLkkb4J/XgYxFtiIj9EdEXEX09PT0plGZmZknSCP/jwIay1+uBZ2ZrI6kFuAh4IYVjm5nZIqQR/o8BV0i6XFIbcD3QX9GmH7ix9Pw9wMMRcU7P38zMlkbVa/tExJSkDwL3AzngixFxWNJeYCAi+oEvAF+WdJRij//6ao9rZmaLl8rCbhFxH3BfxbY9Zc8ngOvSOJaZmVXPV/iamWWQw9/MLIMc/mZmGeTwNzPLIIe/mVkGZfI2jpUOHBli38FBjo2MsWF1J7u3bfZtHs2sqWW+53/gyBB7+g8zNDrBqo5WhkYn2NN/mANHhupdmplZzWQ+/PcdHKQ1JzrbWpCKj605se/gYL1LMzOrmcyH/7GRMTpaczO2dbTmOD4yVqeKzMxqL/Phv2F1J+P5woxt4/kC61d31qkiM7Pay3z47962mXwhGJucIqL4mC8Eu7dtrndpZmY1k/nw376ll707ttLb3c7J8Ty93e3s3bHVs33MrKl5qifFDwCHvZllSeZ7/mZmWeTwNzPLIIe/mVkGOfzNzDLI4W9mlkFVhb+kNZIekPRk6XF1QpvXSvqepMOSfijp16o5ppmZVa/anv+twEMRcQXwUOl1pTHgNyNiK3At8FlJq6o8rpmZVaHa8N8J3FF6fgfwrsoGEfE3EfFk6fkzwBDQU+VxzcysCtWG/yURcQKg9HjeK6UkXQW0AX87y/u7JA1IGhgeHq6yNDMzm82cV/hKehC4NOGtjy/kQJLWAV8GboyIM0ltImI/sB+gr68vFrJ/MzObvznDPyKume09Sc9KWhcRJ0rhnngHFEkXAv8H+EREPLLoas3MLBXVDvv0AzeWnt8I3FvZQFIb8HXgSxHx1SqPZ2ZmKag2/G8D3ibpSeBtpddI6pP0+VKb9wLbgPdL+kHpz2urPK6ZmVVBEY05tN7X1xcDAwP1LsPMbFmRdCgi+uZq5yWdZ3HgyBD7Dg5ybGSMDas72b1ts5d9NrOm4eUdEhw4MsSe/sMMjU6wqqOVodEJ9vQf5sCRxO+zzcyWHYd/gn0HB2nNic62FqTiY2tO7Ds4WO/SzMxS4fBPcGxkjI7W3IxtHa05jo+M1akiM7N0OfwTbFjdyXi+MGPbeL7A+tWddarIzCxdDv8Eu7dtJl8IxianiCg+5gvB7m2b612amVkqHP4Jtm/pZe+OrfR2t3NyPE9vdzt7d2z1bB8zaxqe6jmL7Vt6HfZm1rTc8zczyyCHv5lZBjn8zcwyyOFvZpZBDn8zswxy+JuZZZDD38wsgxz+ZmYZ5PA3M8sgh7+ZWQZVtbyDpDXAXcAm4CngvRExMkvbC4EngK9HxAerOW49+M5eZtZMqu353wo8FBFXAA+VXs/m94H/W+Xx6sJ39jKzZlNt+O8E7ig9vwN4V1IjSb8AXAL8eZXHqwvf2cvMmk214X9JRJwAKD2eMw4iaQXwX4CPVnmsuvGdvcys2cw55i/pQeDShLc+Ps9jfAC4LyKOSZrrWLuAXQAbN26c5+5rb8PqToZGJ+hse/l0+c5eZraczRn+EXHNbO9JelbSuog4IWkdkDQI/gbglyV9AOgC2iSdjohzvh+IiP3AfoC+vr6Y71+i1nZv28ye/sOMTU7R0ZpjPF/wnb3MbFmrdtinH7ix9PxG4N7KBhHx6xGxMSI2AR8BvpQU/I3Md/Yys2ZT7Z28bgPulnQT8DPgOgBJfcAtEXFzlftvGL6zl5k1E0U0zOjKDH19fTEwMFDvMszMlkwa1xNJOhQRfXO18z18zcyWQGWwv2HzGr43+MKM1/c8/jStOc24nmgv1GTUweFvZlYD5WHf1Zbj+X+Y5MKOVlZ1tPLT507zl0+9QG93GxdfsJKh0QluP/C3rLmglYs62gHobGthbHKKfQcHHf5mZo1gob34o8OnmSoEF6wsXig6OjHFCsGp8SnWdrXT2dZC4UxwcizP2q726ePU8noih7+Z2XkkBX15sM+nF184E6wQDI++RHd7K5OFM6wQTBbOTB9nZcsKJqYKM45dy+uJHP5mlmnn68VXDtckBft8evFtuRXkC2emw74tt4LJwhnaci/Ptu9ub2FqLJbseiKH/yJ5lU+z5WehvfjK4ZqkYJ9PL35t10qefnGclhUiIuhub2H49CQXdrQQEYznC7S15Pjt7Rv53uALHB8ZY32Nc8XhvwhnV/lcqm/lzWxxzvel63x68ZXDNXBusM+nF9+SE6s6W+npWsnJ8TyXr+3in1+1JjHoP7RE58bhvwjlq3xC7b+VN7NzVful63x68ZXDNXBusM+3F/8f3nnlOfmwVEGfxOG/CMdGxljV0Tpjm1f5NKutaqdOLqYXXzlckxTsjdCLXwyH/yJ4lU+z2pprbH4xUycX04uvHK45X7A3ctAncfgvglf5NEvXQsfmFzN1crG9+KThmmbg8F+E7Vt62Utx7H8pvpU3ayYL7dWnNXWymXvxi+HwXySv8mk2P9X26iG9qZPN2otfDIe/maWmFr16WB5TJ5cbh7+ZzVu1V8Muple/XKZOLjcOfzNLVIurYRfTq3cvvjYc/inykg+2nC3F1bDV9OotXQ7/lHjJB1tO0hibX8zVsO7VN46qwl/SGuAuYBPwFPDeiBhJaLcR+DywAQjgHRHxVDXHbjRe8sGWi6SOylJdDetefeOotud/K/BQRNwm6dbS648ltPsS8J8i4gFJXcCZhDbLmpd8sEZW3tM/NZ7ngpW5GXeMWuqrYa3+qg3/ncD20vM7gANUhL+kK4GWiHgAICJOV3nMhuQlH6yRnG/8/u9PTjA+WWBlS66qsfksXQ3bjKoN/0si4gRARJyQlPRf/VXAi5K+BlwOPAjcGhGFhLbLlpd8sHpZ6Pj9ypbicE35kI6vhs2eOcNf0oPApQlvfXwBx/hl4HXAzyh+R/B+4AsJx9oF7ALYuHHjPHffGLzkgy2Vaq+YXdu1kmdOjjMxVfDYfIbNGf4Rcc1s70l6VtK6Uq9/HTCU0Ow48P2IGCz9zDeAq0kI/4jYD+wH6Ovri/n9FRqHl3ywtNXiitkLO1p5aarA2GTBY/MZVu2wTz9wI3Bb6fHehDaPAasl9UTEMPAWYKDK45o1paVYB+dsT/+2d7/anZUMqzb8bwPulnQTxSGd6wAk9QG3RMTNEVGQ9BHgIUkCDgF/UuVxzZrC+cJ+Ka6YteyqKvwj4nngrQnbB4Cby14/ALy6mmMtR77i186ncr59Zdj7ilmrJV/hWyO+4tcqVXYGXhybnHFhYGXY+4pZqyWHf434il+ba/z+qefHWL+qfbp9Zdi7V2+15PCvEV/xmy2LmZXTmhPPnnqJCzvagHPD3r16qyWHf434it/mlsasnEu6V3L8xfHzDuG4V2+14vCvEV/x27zm+qJ2vrNyWnIruKKni9UXrPQQji05h3+N+Irf5nK+hdEWOysnXwiHvdWNw7+Gkq749fTPxjfX+H3lwmiea2/LkcN/CXn6Z+Obz1r3lQujeVaOLUcO/yXk6Z+NaaFr3VcujOZZObYcOfyXkKd/1t9Ch3Tg3PH7pIXR3Ku35cbhv4Q8/XPpLXRK5nzWuvfCaNYMVtS7gCzZvW0z+UIxRCKKj57+WTtnx++HRidY1dHKUy+MMTKWp3AmzpmSedbarpUA00M6Y5NTpfH7V9Db3c7J8Ty93e3s3bHVwW/Lmnv+Syhp+ucbNq9h38FBPnHvjzz7JwXVTsn0WveWFQ7/JVY+/dOzf6pTiymZHtKxrPCwTx2Vz/4pX+9l38HBepfW8CqHdM6O3+cLhenzubJlBZR6+lAc0jkTkJM8pGOZ555/HXn2z8J4SqZZehz+dZQ0++e50y8xNlnglz75cKa/A/CUTLPacvjXUeXib8+dfonh05P0drdl7jsAT8k0W1oe86+j7Vt62btj6/R489hkgd7uNtZ2tWfqOwBPyTRbelX1/CWtAe4CNgFPAe+NiJGEdp8C3knxw+YB4HciIqo5drMon/3zS598OPE7gCefPcUN+x9pqsXgPCXTrL6qHfa5FXgoIm6TdGvp9cfKG0h6I/AmXr6B+18AbwYOVHnspjPbdwCjLxVmzGpZbkNBnpJp1niqDf+dwPbS8zsoBvrHKtoE0A60AQJagWerPG5TSroBzMhYnjUXtC6rxeCqHb+f7yqZzfAbkFm9VBv+l0TECYCIOCHpnP8TI+J7kr4DnKAY/n8YEU8k7UzSLmAXwMaNG6ssbflJugL45Hieiy9YOaPdVOEMj/9spCFmBC3m3rWekmlWf5pr6F3Sg8ClCW99HLgjIlaVtR2JiNUVP/9K4L8Cv1ba9ADwsYg4eL7j9vX1xcDAwNx/gyZ3w/5HZgwFnRrPT/eKX9nbNX1HqPe8/jK+N/hCzb8XmK1Xf/Y3leMj46y5oHU62I/8/anir3u5FWzu6QJgcPg0E1MFrlx30fR+h0cnGJsscFFHq3v1ZlWQdCgi+uZqN2fPPyKuOc9BnpW0rtTrXwcMJTT7VeCRiDhd+plvAVcD5w1/K6ocCnp2dAKASy96eUbQcGkoZf3qjumhlY/c81f0dK1k9KWpqj4MzjeEs9h713r83qz+qh326QduBG4rPd6b0OZnwG9J+s8Uh33eDHy2yuNmRuVQUARctqp9ehYMwOjEFFNnzkz/djBVCF4cy3N6YopX9nbN+mFAab/lQzblvz3MNYSz2HvXevzerP7mHPY57w9LFwN3Axsphvx1EfGCpD7gloi4WVIO+BywjeKXv9+OiA/PtW8P+ySrHAYCeOLEKdpy4hW93UBxWCVfOEMAWy69MHGo6OR4HsH0kE35BWYXX7ByXkM4lceB4vDNyFie9as7pnv1J8fz9HSt5PRLUw56sxpLbdjnfCLieeCtCdsHgJtLzwvA7mqOYy9LmhGUWyEu6nz5N4HJwhlEccgFitNFVwgK8fJFU0+/OA4Bl17UARR/e1ghODU+xdqu9nkN4fjetWbLl5d3WGaSZgTtfM0/4p7Hn57+QMitEFOFoKe7OEuo8sMAihdSlf/WN1k4wwoxY2x+riEcz8oxW74c/stQ+VXBZ716/arpD4RNazp5/h8myZVCuvLDACC3QhCaft2WK861L/+AqPxiNins3as3W54c/k2i8gPh7CydpA+D8XyBrpUtCKaDvbu9heHTk1zY0eIhHLMMqOoL31ryF77pKv8wWF8x26f8lpKegWO2vM33C1+Hv5lZE5lv+HtJZzOzDHL4m5llkMPfzCyDHP5mZhnk8Dczy6CGne0jaRj4uyp2sRZ4LqVyloLrrS3XW1uut7YWUu/PRUTPXI0aNvyrJWlgPtOdGoXrrS3XW1uut7ZqUa+HfczMMsjhb2aWQc0c/vvrXcACud7acr215XprK/V6m3bM38zMZtfMPX8zM5uFw9/MLIOWdfhL+qKkIUk/muX97ZJOSvpB6c+epa6xop4Nkr4j6QlJhyX9TkIbSfpvko5K+qGk19ej1lIt86m3Yc6xpHZJfynpr0r1/seENisl3VU6v49K2rT0lU7XMp963y9puOz83lyPWitqykn6vqRvJrzXMOe3rKbz1dtQ51fSU5L+ulTLOcsap5oPEbFs/1C8KfzrgR/N8v524Jv1rrOsnnXA60vPu4G/Aa6saPMO4FuAgKuBRxu83oY5x6Vz1lV63go8Clxd0eYDwB+Xnl8P3NXg9b4f+MN6n9uKmj4MfCXpv3sjnd951ttQ5xd4Clh7nvdTy4dl3fOPiIPAC/WuY74i4kREPF56Pgo8AVxW0Wwn8KUoegRYJWndEpcKzLvehlE6Z6dLL1tLfypnNOwE7ig9vwd4qyRRB/Ost6FIWg+8E/j8LE0a5vzCvOpdblLLh2Ud/vP0htKv1d+StLXexZxV+nX4dRR7e+UuA46VvT5OAwTueeqFBjrHpV/xfwAMAQ9ExKznNyKmgJPAxUtb5cvmUS/APyv9in+PpA1LXGKlzwL/Djgzy/sNdX6Zu15orPMbwJ9LOiRpV8L7qeVDs4f/4xTXuXgN8N+Bb9S5HgAkdQH/C/jXEXGq8u2EH6lrb3COehvqHEdEISJeC6wHrpL08xVNGur8zqPe/w1siohXAw/ycq96yUn6p8BQRBw6X7OEbXU5v/Ost2HOb8mbIuL1wNuB35a0reL91M5vU4d/RJw6+2t1RNwHtEpaW8+aJLVSDNL/GRFfS2hyHCjvfawHnlmK2pLMVW8jnuNSLS8CB4BrK96aPr+SWoCLaIChw9nqjYjnI+Kl0ss/AX5hiUsr9yZgh6SngDuBt0j6HxVtGun8zllvg51fIuKZ0uMQ8HXgqoomqeVDU4e/pEvPjjdKuori3/f5OtYj4AvAExHxB7M06wd+s/St/tXAyYg4sWRFlplPvY10jiX1SFpVet4BXAMcqWjWD9xYev4e4OEofZO21OZTb8V47g6K37vURUT8bkSsj4hNFL/MfTgi3lfRrGHO73zqbaTzK+kCSd1nnwP/BKicyZhaPrRUVW2dSfozirNN1ko6DvwexS/NiIg/pviP719JmgLGgevr9Q+x5E3AbwB/XRrnBfj3wEaYrvk+it/oHwXGgH9RhzrPmk+9jXSO1wF3SMpR/BC6OyK+KWkvMBAR/RQ/zL4s6SjFHun1daoV5lfvhyTtAKYo1vv+ulU7iwY+v4ka+PxeAny91JdqAb4SEd+WdAuknw9e3sHMLIOaetjHzMySOfzNzDLI4W9mlkEOfzOzDHL4m5mlSHMsOLmI/X1b0ouVC9NJury0eN6TKi6m17aQ/Tr8zczS9aece3FhNT5Nccp1pU8Cn4mIK4AR4KaF7NThb2aWoqQFJyW9otSDPyTpu5K2LGB/DwGjFfsT8BaKi+dBcVmKdy2kzmV9kZeZ2TKxH7glIp6U9IvA5yiG92JdDLxYWjwPFrHAm8PfzKyGSgsjvhH4atnq1itL770b2JvwY09HxK+cb7cJ2xZ0xa7D38ystlZQ7KW/tvKN0mKJSQs8zuU5imv5t5R6/wte4M1j/mZmNVRaBv2nkq6D6VsxvqbKfQbwHYpra0FxMb17F7IPr+1jZpai8gUngWcpLjj5MPBHFBfzawXujIik4Z6k/X0X2AJ0UVwx96aIuF/SZopLVa8Bvg+8r2x56rn36/A3M8seD/uYmWWQw9/MLIMc/mZmGeTwNzPLIIe/mVkGOfzNzDLI4W9mlkH/H3XWZzsJSeSiAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "R=[i*Ao/20 for i in range(30,100)]\n",
    "Y=[Coul_Buck(r, (r,0,0), Na,Cl,NaCl)[0] for r in R]\n",
    "sns.regplot(x=R,y=Y,fit_reg=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x120856eb8>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEQCAYAAABbfbiFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGzZJREFUeJzt3X+QXeV93/H3x6sV2jWytVi7jrBEQAmtEjEF3A0B41EFoTGRXfGjpgMTPKQhIzl1WjceO4Smo8RyMhPqpnjaxM4qNjFpbMChVlEBg/m1s7gewCtbCMkSgSzEksV4ZbNg0ZVgtfr2j3MW7u65u/es7r177o/Pa+bO3vuc557z1dG993ufH+e5igjMzMxKva3oAMzMrPE4OZiZWYaTg5mZZTg5mJlZhpODmZllODmYmVlG0yYHSbdJGpW0p0b7e0DSK5LunVH+K5K+K2mXpG9J+vlaHM/MrJE1bXIAvgxcXsP9fRb4SJnyLwC/HhHnAV8F/nMNj2lm1pCaNjlExBDwcmmZpJ9LWwA7JT0uac089vcIcKTcJuAd6f13AodONmYzs2axqOgAamwb8NGIeE7SLwOfBy6tcp+/Bdwv6SjwU+DCKvdnZtbwWiY5SDoVeB/wd5Kmik9Jt10NbC3ztB9GxAcq7Pp3gQ0R8aSkTwH/jSRhmJm1rJZJDiRdZK+kYwPTRMTXga/Pd4eSeoFzI+LJtOgu4IGqojQzawJNO+YwU0T8FHhB0jUASpxb5W7HgHdK+ifp438J7Ktyn2ZmDU/NuiqrpDuA9cBy4EfAHwKPkswuWgF0AndGRLnupHL7exxYA5wK/AS4MSIelHQVSZfUCZJk8ZsRMVLbf42ZWWNp2uRgZmb10zLdSmZmVjtNOSC9fPnyOPPMM4sOw8ysqezcufPHEdGbp25TJoczzzyT4eHhosMwM2sqkv4xb113K5mZWYaTg5mZZTg5mJlZhpODmZllVJ0cJF0jaa+kE5L6Z6mzRNJTkp5O6366ZNtZkp6U9JykuyQtrjYmMzOrTi1aDnuAq4GhOeq8DlwaEecC5wGXS5pa3fQW4NaIOJvkCuQbaxCTmZlVoerkEBH7IuLZCnUiIl5LH3amt1CyfOqlwN3pttuBK6uNyczMqrNg1zlI6gB2Aj8P/EW6BPZykpVUj6fVDgLvmeX5m4BNAGecccYCRGxm1hgG948yMDTCgbFxVvV0s3ndatav6avrMXO1HCQ9LGlPmdsVeQ8UEZPpctorgQsknQOoXNVZnr8tIvojor+3N9cFfmZmTW9w/yhbduxl9MgxlnV1MnrkGFt27GVw/2hdj5ur5RARl9XqgBHxiqRBkt9//jNgmaRFaethJf4ZTjOzNw0MjdDZIboXJx/X3YsXMf7GcQaGRuraeliQqaySeiUtS+93AZcB+yNZEvYx4MNp1RuAexYiJjOzZnBgbJyuzo5pZV2dHRwcG6/rcWsxlfUqSQeBi4D7JD2Ylp8u6f602grgMUm7ge8AD0XEvem2m4BPSHoeeBfwpWpjMjNrFat6ujk6MTmt7OjEJCt7uut63KoHpCNiO7C9TPkhYEN6fzdw/izPHwEuqDYOM7NWtHndarbs2Mv4G8fp6uzg6MQkE5PB5nWr63pcXyFtZtbA1q/pY+vGtfQtXcKrRyfoW7qErRvX1n22UlMu2W1m1k7Wr+mrezKYyS0HMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzs4yqkoOkayTtlXRCUv8sdZZIekrS02ndT5ds+7KkFyTtSm/nVROPmZnVxqIqn78HuBoYmKPO68ClEfGapE7gW5K+ERFPpNs/FRF3VxmHmZnVUFXJISL2AUiaq04Ar6UPO9NbVHNcMzOrrwUZc5DUIWkXMAo8FBFPlmz+E0m7Jd0q6ZQ59rFJ0rCk4cOHD9c9ZjOzdlYxOUh6WNKeMrcr8h4kIiYj4jxgJXCBpHPSTTcDa4BfAk4DbppjH9sioj8i+nt7e/Me2szMTkLFbqWIuKxWB4uIVyQNApcDeyLipXTT65L+GvhkrY5lZtYMBvePMjA0woGxcVb1dLN53WrWr+krOqz6dytJ6pW0LL3fBVwG7E8fr0j/CriSZIDbzKwtDO4fZcuOvYweOcayrk5Gjxxjy469DO4fLTq0qqeyXiXpIHARcJ+kB9Py0yXdn1ZbATwmaTfwHZIxh3vTbV+R9AzwDLAc+ONq4jEzayYDQyN0dojuxYuQkr+dHWJgaKTo0KqerbQd2F6m/BCwIb2/Gzh/ludfWs3xzcya2YGxcZZ1dU4r6+rs4ODYeEERvcVXSJuZFWRVTzdHJyanlR2dmGRlT3dBEb3FycHMrCCb161mYjIYf+M4Ecnficlg87rVRYfm5GBmVpT1a/rYunEtfUuX8OrRCfqWLmHrxrUNMVup2uUzmkqjThkzs/a1fk1fQ34OtU3LoZGnjJmZNZq2SQ6NPGXMzKzRtE1yODA2Tldnx7SyRpkyZmbWaNomOTTylDEzs0bTNsmhkaeMmZk1mrZJDo08ZczMrNG01VTWRp0yZmbWaNqm5WBmZvk5OZiZWYaTg5mZZTg5mJlZhpODmZllODmYmVmGk4OZmWU4OZiZWYaTg5mZZTg5mJlZhpODmZllODmYmVmGk4OZmWU4OZiZWYaTg5mZZTg5mJlZhpODmZllODmYmVlGW/1MqJnZQhrcP8rA0AgHxsZZ1dPN5nWrm+aniqtuOUi6RtJeSSck9Veo2yHpe5LuLSk7S9KTkp6TdJekxdXGZGZWtMH9o2zZsZfRI8dY1tXJ6JFjbNmxl8H9o0WHlkstupX2AFcDQznqfhzYN6PsFuDWiDgbGANurEFMZmaFGhgaobNDdC9ehJT87ewQA0MjRYeWS9XJISL2RcSzlepJWgl8EPhiSZmAS4G706LbgSurjcnMrGgHxsbp6uyYVtbV2cHBsfGCIpqfhRyQ/hzwe8CJkrJ3Aa9ExPH08UHgPeWeLGmTpGFJw4cPH65vpGZmVVrV083RiclpZUcnJlnZ011QRPOTKzlIeljSnjK3K3I+/0PAaETsnLmpTPUot4+I2BYR/RHR39vbm+ewZmaF2bxuNROTwfgbx4lI/k5MBpvXrS46tFxyzVaKiMuqPM7FwEZJG4AlwDsk/S3wEWCZpEVp62ElcKjKY1WlmWcXmFnjWL+mj60kYw8Hx8ZZ2WSfJwsylTUibgZuBpC0HvhkRFyfPn4M+DBwJ3ADcM9CxFTO1OyCzg5Nm12wFZrmP9TMGsf6NX1N+9lRi6msV0k6CFwE3CfpwbT8dEn359jFTcAnJD1PMgbxpWpjOlnNPrvAzKxWqm45RMR2YHuZ8kPAhjLlg8BgyeMR4IJq46iFA2PjLOvqnFbWTLMLzMxqxctnlGj22QVmZrXi5FCi2WcXmJnVipNDifVr+ti6cS19S5fw6tEJ+pYuYevGtU07oGRmdrK88N4MzTy7wMysVtxyMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswwnBzMzy3ByMDOzDK/KamZ2Egb3jzIwNMKBsXFW9XSzed3qllrR2S0HM7N5Gtw/ypYdexk9coxlXZ2MHjnGlh17Gdw/WnRoNePkYGY2TwNDI3R2iO7Fi5CSv50dYmBopOjQasbJwcxsng6MjdPV2TGtrKuzg4Nj4wVFVHtODmZm87Sqp5ujE5PTyo5OTLKyp7ugiGrPyeEkDO4f5bptT/D+Wx7lum1PtFQ/o5lVtnndaiYmg/E3jhOR/J2YDDavW110aDXj5DBP7TAQZWZzW7+mj60b19K3dAmvHp2gb+kStm5c21KzlTyVdZ5KB6IAuhcvYvyN4wwMjbTUC8PM5rZ+TV9Lv+fdcpindhiIMjNzcpindhiIMjNzcpindhiIMjNzcpindhiIMjPzgPRJaPWBKDOzqloOkq6RtFfSCUn9Fep2SPqepHtLyr4s6QVJu9LbedXEY2ZmtVFty2EPcDUwkKPux4F9wDtmlH8qIu6uMg4zM6uhqloOEbEvIp6tVE/SSuCDwBerOZ6ZmS2MhRqQ/hzwe8CJMtv+RNJuSbdKOmW2HUjaJGlY0vDhw4frFqiZmeVIDpIelrSnzO2KPAeQ9CFgNCJ2ltl8M7AG+CXgNOCm2fYTEdsioj8i+nt7e/Mc2szMTlLFMYeIuKzKY1wMbJS0AVgCvEPS30bE9RHxUlrndUl/DXyyymOZmVkN1L1bKSJujoiVEXEmcC3waERcDyBpRfpXwJUkA9xmZlawaqeyXiXpIHARcJ+kB9Py0yXdn2MXX5H0DPAMsBz442riMTOz2lBEFB3DvPX398fw8HDRYZhZixrcP8rA0AgHxsZZ1dPN5nWrW+LCV0k7I2LOa9KmePkMM7MS/s2WhJODmVmJ0t9skZK/nR1iYGik6NAWlJODmVkJ/2ZLwgvv1UGr9leatYNVPd2MHjn25q89Qnv+ZotbDjXm/kqz5ubfbEk4OdSY+yvNmpt/syXhbqUaOzA2zrKuzmll7dhfadbM/JstbjnUnH9j2sxagZNDjbm/0sxagZNDjbm/0sxagccc6sD9lWbW7NxyMDOzDCcHMzPLcHIwM7MMJwczM8twcjAzswzPVjKztuKFMfNxy8HM2oYXxszPycHM2oYXxszP3UoFcdPWbOF5Ycz83HIogJu2ZsXwwpj5OTkUwE1bs2J4Ycz8nBwK4N+oNSuGF8bMz2MOBfBv1JoVxwtj5uOWQwHctDWzRufkUAA3bc2s0blbqSBu2ppZI3PLwczMMpwczMwsw8nBzMwyqhpzkHQN8EfALwAXRMTwLPVeBI4Ak8DxiOhPy08D7gLOBF4E/k1EjFUTk5m1Ly9LUzvVthz2AFcDQznqXhIR500lhtTvA49ExNnAI+ljM7N587I0tVVVyyEi9gFIOtldXAGsT+/fDgwCN1UTUyvxtyCz/EqXpQHoXryI8TeOMzA04vfNSVioMYcAvilpp6RNJeXvjoiXANK/s/4PStokaVjS8OHDh+scbvH8LchsfrwsTW1VTA6SHpa0p8ztinkc5+KIeC/wa8DHJK2bb6ARsS0i+iOiv7e3d75PbzpenM9sfrziam1VTA4RcVlEnFPmdk/eg0TEofTvKLAduCDd9CNJKwDSv/5anPK3ILP58bI0tVX3biVJb5e0dOo+8KskA9kAO4Ab0vs3ALkTTqvztyCz+fGyNLVV7VTWq4D/AfQC90naFREfkHQ68MWI2AC8G9ieDlovAr4aEQ+ku/hT4GuSbgR+AFxTTTytZPO61WzZsZfxN47T1dnB0YlJfwsyq8DL0tSOIqLoGOatv78/hofLXlLRUqZmKx0cG2elZyuZWZUk7ZxxOcGsvPBeA/O3IDMripODmTUNX/uzcLy2kpk1BV/7s7Dccmhy/iZl7cJXQC8stxyamL9JWTvxtT8Ly8mhifkqamsnvvZnYTk5NDF/k7J24iugF5aTQxPzNylrJ74CemF5QLqJ+Spqaze+9mfhuOXQxPxNyszqxS2HJudvUtYqPC27sTg5tAG/6azRTU3L7uzQtGnZW8Gv1YK4W6nF+VoIawaelt14nBxanN901gw8LbvxODm0OL/prBl4WnbjcXJocX7TWTPwBW6Nx8mhxflNZ41gcP8o1217gvff8ijXbXsiM+bladmNx78E1wby/KKcZzRZvZTORCq9WNMf/gvPvwRn01S6FsLTCK2evNR2c3K3knlGk9WVJ0U0JycH85vX6sqTIpqTu5WMVT3djB459mazH8q/eT0uYeVUel14gcjm5JaD5ZrR5CutrZw8rwvPRGpObjlY8uaFOWc0eVDRysn7uvACkc3HycGAym/eA2PjLOvqnFbmcYnWV6nLyK+L1uXkYLl4XKL95JninPd1Yc3HYw6Wi8cl2k+eKc6+Ar91ueVgudRyXMKti8ZQiy6jPK8La05ODpZbLcYlfDV2Y6hll5EHm1uTu5WsZvJc7JT3auxKC7XZ3CqdP3cZWSVVJQdJ10jaK+mEpFkXc5L0oqRnJO2SNFxS/keSfpiW75K0oZp4rFh5PkzyXI3tsYu5Vfrgz3P+8vw/+PqE9lZtt9Ie4GpgIEfdSyLix2XKb42I/1plHNYA8vQ/5+mqqOXYRbONb1SKN093UJ7z5y4jq6Sq5BAR+wAk1SYaa3qVPkzyLKVQq7GLvOMbeT6Qa5GEFuqDP8/585IWVslCjTkE8E1JOyVtmrHtdyTtlnSbpJ7ZdiBpk6RhScOHDx+ub7RWN3m6Kmo1dpGnTqUumDxdNLWqkyfePN1Bec6fu4yskootB0kPAz9TZtMfRMQ9OY9zcUQcktQHPCRpf0QMAV8APkOSPD4D/Bnwm+V2EBHbgG2Q/NhPzuNaA1qo1kWeOpW+ief5pl6rOnnizdMdlLdV4C4jm0vFlkNEXBYR55S55U0MRMSh9O8osB24IH38o4iYjIgTwF9NlVt7q1XrIk+dSt/E83xTr1WdPPHmGfR3q8Bqoe7XOUh6O/C2iDiS3v9VYGu6bUVEvJRWvYpkgNusJq2LPHUqfRPP8029VnXyxJv3ojO3Cqxa1U5lvUrSQeAi4D5JD6blp0u6P632buBbkp4GngLui4gH0m3/JZ3iuhu4BPjdauKx9pHn23GeOpW+ief5pl6rOnm/8a9f08cdmy7k8Zsu5Y5NFzoJWF0oovm67/v7+2N4eLhyRbMcpmYRzfZNvNL2WtYxqydJOyNi1mvSptV1cjAzaw/zSQ5ePsPMzDKcHMzMLMPJwczMMpwczMwsw8nBzMwynBzMzCzDycHMzDKa8joHSYeBf6xiF8uBcr8t0aiaLV5wzAuh2eKF5ou52eKFuWP+2YjozbOTpkwO1ZI0nPdCkEbQbPGCY14IzRYvNF/MzRYv1C5mdyuZmVmGk4OZmWW0a3LYVnQA89Rs8YJjXgjNFi80X8zNFi/UKOa2HHMwM7O5tWvLwczM5uDkYGZmGS2THCStkvSYpH2S9kr6eJk6kvTfJT0vabek95Zsu0HSc+nthgaK+dfTWHdL+rakc0u2vZj+kt4uSQvyAxc5Y14v6dU0rl2StpRsu1zSs+n/we83SLyfKol1j6RJSael24o4x0skPSXp6TTmT5epc4qku9Lz+KSkM0u23ZyWPyvpAw0S7yckfT99HT8i6WdLtk2WnP8d9Y53HjH/hqTDJbH9Vsm2Ij4v8sR8a0m8fy/plZJt8zvPEdESN2AF8N70/lLg74FfnFFnA/ANQMCFwJNp+WnASPq3J73f0yAxv28qFuDXpmJOH78ILG/A87weuLfMczuAfwBWA4uBp2c+t4h4Z9T/V8CjBZ9jAaem9zuBJ4ELZ9T5d8BfpvevBe5K7/9iel5PAc5Kz3dHA8R7CdCd3v/tqXjTx68t5PmdR8y/Afx5mecW9XlRMeYZ9f89cNvJnueWaTlExEsR8d30/hFgH/CeGdWuAP4mEk8AyyStAD4APBQRL0fEGPAQcHkjxBwR305jAngCWFnvuOaS8zzP5gLg+YgYiYg3gDtJ/k/q5iTivQ64o54xVZK+Pl9LH3amt5kzR64Abk/v3w38iiSl5XdGxOsR8QLwPMl5LzTeiHgsIsbTh43wOs5zjmdT1OfFfGOu6rXcMsmhVNrEPp8ks5Z6D3Cg5PHBtGy28gUzR8ylbiRp+UwJ4JuSdkraVL/oyqsQ80Vp8/cbktamZYWe50rnWFI3yZv8f5UUF3KOJXVI2gWMknwQzfpajojjwKvAuyjoHOeIt9TM1/ESScOSnpB0ZV0DLZEz5n+ddoXdLWlVWlbY6zjveU677c4CHi0pntd5XlSTiBuIpFNJ3tz/MSJ+OnNzmafEHOULokLMU3UuIXlTvb+k+OKIOCSpD3hI0v6IGKp/xBVj/i7JGi6vSdoA/G/gbAo8z3nOMUmX0v+NiJdLygo5xxExCZwnaRmwXdI5EbGnpEpDvZZzxAuApOuBfuBflBSfkZ7j1cCjkp6JiH9ogJj/D3BHRLwu6aMkLbVLKfB1nPc8k3Q13p3WnzKv89xSLQdJnSQfAF+JiK+XqXIQWFXyeCVwaI7yussRM5L+GfBF4IqI+MlUeUQcSv+OAtupc/dBSTxzxhwRP51q/kbE/UCnpOUUdJ7znOPUtcxohhd1jkuO/wowSLbb4s1zKWkR8E7gZQp8LcOc8SLpMuAPgI0R8XrJc6bO8Uj63PMXItaS45eNOSJ+UhLnXwH/PL1f6DlOY5v1PKfmei3nO8/zGaBo5BtJNv8b4HNz1Pkg0wekn4q3BpheIBlc6knvn9YgMZ9B0m/8vhnlbweWltz/NnB5g8T8M7x1geUFwA/S5y0iGbw7i7cGpNcWHW9ab+rD9e0NcI57gWXp/S7gceBDM+p8jOkD0l9L769l+oD0CPUfkM4T7/kkg+NnzyjvAU5J7y8HnqPOkxTmEfOKkvtXAU+k94v6vKgYc7rtn5JMpFA157mVupUuBj4CPJP2yQH8J5IPVyLiL4H7SWYsPQ+MA/823faypM8A30mftzWmdy0UGfMWkr7kzyfjjRyPZMXFd5M0KyH50P1qRDzQIDF/GPhtSceBo8C1kbwqj0v6HeBBkplLt0XE3gaIF5I3/zcj4v+VPLeoc7wCuF1SB0nr/msRca+krcBwROwAvgT8T0nPkyS1a9N/z15JXwO+DxwHPhbTuxaKivezwKnA36Xn8wcRsRH4BWBA0on0uX8aEd+vc7x5Y/4PkjaSnMeXSWYvFfl5kSdmSAai70zfc1PmfZ69fIaZmWW01JiDmZnVhpODmZllODmYmVmGk4OZmWU4OZiZLQBJt0kalVTuorWT2d8Dkl6RdO+M8rOULMb4nJLFGRefzP6dHMzMFsaXqe0aTJ8lmaY90y3ArRFxNjBGsrLCvDk5mJktgEiWXZl2PYSkn0tbADslPS5pzTz29whwZMb+RLLEx91p0e3ASa1X1UoXwZmZNZttwEcj4jlJvwx8nuTD/WS9C3glksUYoYpFAZ0czMwKkC4G+T7eumockmVPkHQ1sLXM034YEXP9gFPNFgV0cjAzK8bbSL7lnzdzQyQLRM61SORsfkzyOzWL0tbDSS8K6DEHM7MCRLJ0/AuSroE3f8b43ApPq7TPAB4jWd8M4AbgnpPZl9dWMjNbAJLuIPkJ3eXAj4A/JPkxni+QLKrXSbJgXrnupHL7exxYQ7Kg4U+AGyPiwfT3Gu4kWT32e8D1UbJEeu54nRzMzGwmdyuZmVmGk4OZmWU4OZiZWYaTg5mZZTg5mJlZhpODmZllODmYmVnG/weHdH2wSs0SHAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "R=[i*Ao/20 for i in range(40,75)]\n",
    "Y=[Coul_Buck(r, (r,0,0), Na,Cl,NaCl)[0]+Coul_Buck(Na.latt_const-r, (Na.latt_const-r,0,0), Na,Cl,NaCl)[0] for r in R]\n",
    "sns.regplot(x=R,y=Y,fit_reg=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x1213f4eb8>"
      ]
     },
     "execution_count": 74,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEQCAYAAABLMTQcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAHBxJREFUeJzt3X+Q3Hd93/Hn+3b3fulOlqy7AyOdIisoVWVqMNwYEzpCCUmQCSOlLiQWNSmMiUSLoZ2macy0FRkxmZowHUIzBqy6jgNt7ICHpCojcCY4HnkaBD4BoZExjkZ20FnYd5JO0p3v5+6++8funb/31f747t13b3e/93rMaHS7348++9F3dl766PP9fN9fc3dERCRZ2ho9ABERiZ/CXUQkgRTuIiIJpHAXEUkghbuISAIp3EVEEqih4W5mD5nZqJn9XUz9fdPMLpvZ10PvP2xmz5vZD4q/3hTH54mINKtGz9wfBvbG2N9ngA+UOfY77v6m4q8fxPiZIiJNp6Hh7u4ngEvB98zsZ4sz8FNm9pSZ7ayhv28BE3GPU0Sk1TR65l7KUeBj7v4W4N8Dn4+p3983sx+a2WfNrCOmPkVEmlK60QMIMrMe4OeBr5rZwtsdxWN3AEdK/LEX3f1dVbr+BPAS0E7hH4/fLdOXiEgiNFW4U/ifxGV3v+aCp7t/Dfjacjp1958Wf5w1sz+m8D8CEZHEaqplGXe/CjxvZu8DsII3rrRfM7thoT/g14BYdueIiDQra2RVSDN7BNgD9AEvA58EngC+ANwAZIBH3T3SEoqZPQXsBHqAi8Dd7v64mT0B9AMG/AD4iLtPxvu3ERFpHg0NdxERqY+mWpYREZF4NOyCal9fn2/btq1RHy8i0pJOnTp1wd37q7VrWLhv27aN4eHhRn28iEhLMrN/iNJOyzIiIgmkcBcRSaCq4V6tcqOZ/Yvibf0/NLO/iWNfuoiIrEyUmfvDVK7c+DzwDne/GfgUhdv7RUSkgapeUHX3E2a2rcLxvwm8PAlsWfmwRERkJeLeLXM38I1yB83sIHAQYOvWrTF/tIhI83ry2VEeOHGWc+NTDG7s5tDu7ezZOVC3z4vtgqqZ/QKFcP/dcm3c/ai7D7n7UH9/1W2aIiKJ8OSzoxw+dprRiRk2dGUYnZjh8LHTPPnsaN0+M5ZwN7ObgQeB/e5+MY4+RUSS4oETZ8mkjO72NGaF3zMp44ETZ+v2mSsOdzPbSqEU7wfc/bmVD0lEJFnOjU/RlUktea8rk2JkfKpun1l1zT1YudHMRihUbswAuPsXgcPAJuDzxQdsZN19qF4DFhFpNYMbuxmdmKG7/dXInZ7PsWVjd90+M8pumQNVjn8Y+HBsIxIRSZhDu7dz+NhppuaydGVSTM/nmM85h3Zvr9tn6g5VEZE627NzgCP7bmKgt5Mr0/MM9HZyZN9Ndd0t02yP2RMRSaQ9OwfqGuZhmrmLiCSQwl1EJIEU7iIiCaRwFxFJIIW7iEgCKdxFRBJI4S4ikkAKdxGRBNJNTCIidbLaNdyDNHMXEamDRtRwD1K4i4jUQSNquAcp3EVE6qARNdyDFO4iInUwuLGb6fnckvfqXcM9SOEuIlIHh3ZvZz7nTM1lcS/8Xu8a7kEKdxGROmhEDfcgbYUUEamT1a7hHqSZu4hIAincRUQSSOEuIpJACncRkQRSuIuIJJDCXUQkgRTuIiIJVHWfu5k9BLwHGHX3N5Q4bsDngHcDU8AH3f17cQ9URKTZNbLEb1iUmfvDwN4Kx28HdhR/HQS+sPJhiYi0lkaX+A2rGu7ufgK4VKHJfuBLXnAS2GBmN8Q1QBGRVtDoEr9hcay5bwbOBV6PFN+7hpkdNLNhMxseGxuL4aNFRJpDo0v8hsUR7lbiPS/V0N2PuvuQuw/19/fH8NEiIs2h0SV+w+II9xFgMPB6C3A+hn5FRFpGo0v8hsUR7seA37SC24Ar7v7TGPoVEWkZjS7xGxZlK+QjwB6gz8xGgE8CGQB3/yJwnMI2yDMUtkJ+qF6DFRFpZo0s8RtWNdzd/UCV4w58NLYRiYjIiukOVRGRBFK4i4gkkMJdRCSBFO4iIgmkcBcRSaCqu2VERKS8ZqoEGaSZu4jIMjVbJcgghbuIyDI1WyXIIIW7iMgyNVslyCCFu4jIMjVbJcgghbuIyDI1WyXIIIW7iMgyNVslyCBthRQRWYFmqgQZpJm7iEgCKdxFRBJI4S4ikkAKdxGRBNIFVRGRGjVrPZkgzdxFRGrQzPVkghTuIiI1aOZ6MkEKdxGRGjRzPZkghbuISA2auZ5MkMJdRKQGzVxPJkjhLiJSg2auJxMUaSukme0FPgekgAfd/b7Q8a3AnwAbim3udffjMY9VRKQpNGs9maCqM3czSwH3A7cDu4ADZrYr1Ow/AV9x91uAO4HPxz1QERGJLsqyzK3AGXc/6+5zwKPA/lAbB9YXf74OOB/fEEVEpFZRwn0zcC7weqT4XtDvAXeZ2QhwHPhYqY7M7KCZDZvZ8NjY2DKGKyIiUURZc7cS73no9QHgYXf/r2b2NuDLZvYGd88v+UPuR4GjAENDQ+E+RESaUiuUGwiLMnMfAQYDr7dw7bLL3cBXANz920An0BfHAEVEGqlVyg2ERQn3p4EdZnajmbVTuGB6LNTmJ8A7AczsH1MId627iEjLa5VyA2FVw93ds8A9wOPAjyjsijltZkfMbF+x2W8Dv2Vmfws8AnzQ3bXsIiItr1XKDYRF2ude3LN+PPTe4cDPzwBvj3doIiKNN7ixm9GJGbrbX43LZiw3EKY7VEVEKmiVcgNhCncRkQpapdxAmJ7EJCJSRSuUGwjTzF1EJIEU7iIiCaRwFxFJIK25i4iU0IolB4I0cxcRCWnVkgNBCncRkZBWLTkQpHAXEQlp1ZIDQQp3EZGQwY3dTM/nlrzXCiUHghTuIiIhrVpyIEjhLiIS0qolB4K0FVJEpIRWLDkQpJm7iEgCaeYuIkLr37QUppm7iKx5SbhpKUzhLiJrXhJuWgpTuIvImpeEm5bCFO4isuYl4aalMIW7iKx5SbhpKUzhLiJrXhJuWgrTVkgRaUlxb11s9ZuWwhTuItJyFrYuZlK2ZOviEagpoJO2tz0o0rKMme01sx+b2Rkzu7dMm183s2fM7LSZ/Wm8wxQReVUcWxeTuLc9qGq4m1kKuB+4HdgFHDCzXaE2O4BPAG9395uAf1uHsYqIAPFsXUzi3vagKDP3W4Ez7n7W3eeAR4H9oTa/Bdzv7uMA7p6Mf/pEpCnFsXUxiXvbg6KE+2bgXOD1SPG9oJ8Dfs7M/q+ZnTSzvXENUEQkLI6ti0nc2x4UJdytxHseep0GdgB7gAPAg2a24ZqOzA6a2bCZDY+NjdU6VhERIJ6ti0nc2x4UZbfMCDAYeL0FOF+izUl3nweeN7MfUwj7p4ON3P0ocBRgaGgo/A+EiEhkK926uGfnAEcorL2PjE+xJWG7ZaKE+9PADjO7EXgRuBN4f6jNX1CYsT9sZn0UlmmScVVCRBIraXvbg6qGu7tnzewe4HEgBTzk7qfN7Agw7O7Hisd+xcyeAXLA77j7xXoOXESkVkne1x5m7o1ZHRkaGvLh4eGGfLaIrD3BG5+6Mimm53PM57zlygyY2Sl3H6rWTrVlRGRNSPq+9jCFu4isCUnf1x6mcBeRNSHp+9rDFO4isiYkfV97mKpCikiiBXfI9LSnMDOuTM8nbl97mMJdRBIrXBq4sEMmz6f2vyGxob5AyzIiklhrbYdMkMJdRBJrre2QCVK4i0hirbUdMkFacxeRxFm4iPrcy1eZnM1x/boMm9Z1LN6VmtQdMkEKdxFJlOBF1Buu6+LC5CyXXplnPufsGOhN9A6ZIIW7iLS84HbHq9PzrOtIcV1XJwD9vZ2s60gz0NvJIwdva/BIV4/CXURaWni740tXZpiey9GRTtHbmQHWzkXUIF1QFZGWFt7u2JFuA4OxidnFNmvlImqQwl1EWlp4u2NfTwcAM9ncmigzUI6WZUSkpQ1u7GZ0Yobu9kKcre/KMJvNMTWXWxNlBspRuItISzu0ezuHj51mai67+BCO9nSK++64ec0FepCWZUSkpe3ZOcCRfTcx0NvJlel5Bno7W+7pSvWgmbuItLwkP+h6uTRzFxFJIIW7iEgCKdxFRBJI4S4ikkAKdxGRBFK4i4gkUKRwN7O9ZvZjMztjZvdWaPdeM3MzG4pviCIiUquq4W5mKeB+4HZgF3DAzHaVaNcLfBz4TtyDFBGR2kSZud8KnHH3s+4+BzwK7C/R7lPAHwAzMY5PRESWIUq4bwbOBV6PFN9bZGa3AIPu/vVKHZnZQTMbNrPhsbGxmgcrIiLRRAl3K/GeLx40awM+C/x2tY7c/ai7D7n7UH9/f/RRiohITaKE+wgwGHi9BTgfeN0LvAF40sxeAG4DjumiqohI40QJ96eBHWZ2o5m1A3cCxxYOuvsVd+9z923uvg04Cexz9+G6jFhERKqqGu7ungXuAR4HfgR8xd1Pm9kRM9tX7wGKiEjtIpX8dffjwPHQe4fLtN2z8mGJiMhK6A5VEZEEUriLiCSQwl1EJIEU7iIiCaRwFxFJIIW7iEgCKdxFRBJI4S4ikkAKdxGRBFK4i4gkkMJdRCSBFO4iIgmkcBcRSaBIVSFFRBrhyWdHeeDEWc6NTzG4sZtDu7ezZ+dAo4fVEhTuItJUFgL9uZevMjmb4/p1GTat62B0YobDx05zBBTwEWhZRkSaxpPPjnL42GlGJ2aYmc+Td+fi5DyTs1m629NkUsYDJ842epgtQTN3EWmo4NLL1el51nWkuK6rk7lcnpQZDoxNzNLbmaErk2JkfKrRQ24JCncRaZiFmXomZWzoyvDSlRmm53J0pFO0p9rI5hxrg7lcHoDp+RxbNnY3eNStQcsyItIwD5w4SyZldLenMTM60m1ghZl6X08HeZxc3sm0GVNzWeZzzqHd2xs97JagmbuIrLqFpZjvvnCJjpQxsL6T3s4MfT0dnL8yzUw2R29nmk3Zdsan5unuSDPQ26ndMjVQuIvIqgouxXSm25jL5Tl/eYbXbYD1XRlmszmm5nJcmZ7nxr4e7lOgL4vCXURWVXApZmGm7jijV2dItRnt6RT33XGzAn2FFO4iUnfBHTFjE7O8dn0HUJipA1yYnGUmm9fSS4wU7iJSV+EdMRcmZ3nx8gxmRm9nhvVdGdIpY6C3k0cO3tbo4SZGpN0yZrbXzH5sZmfM7N4Sx/+dmT1jZj80s2+Z2c/EP1QRaUXhHTGv6e0E4KUrM7h75F0wTz47yoGjJ/mnn36CA0dP8uSzo6sx/JZVNdzNLAXcD9wO7AIOmNmuULPvA0PufjPwGPAHcQ9URFrLQhh/94VL/PTyNBMz80BhKWbzhk4cuDI9z0BvJ0f23VRxKSZ45+qGrsxiKQIFfHlRlmVuBc64+1kAM3sU2A88s9DA3f860P4kcFecgxSR1lJpR0xvZ4Z0qo03b90YeRkmOPsH6G5PMzWX5YETZ7U+X0aUZZnNwLnA65Hie+XcDXyj1AEzO2hmw2Y2PDY2Fn2UItJSwjtigMUdMcu5Genc+BRdmdSS91SKoLIoM3cr8Z6XbGh2FzAEvKPUcXc/ChwFGBoaKtmHiLSuUjcnxbEjZnBjN6MTM4szd1ApgmqihPsIMBh4vQU4H25kZr8E/EfgHe4+G8/wRKRVVLs5aSU7Yg7t3s7hY6eZmsvSlUkxPZ9TKYIqooT708AOM7sReBG4E3h/sIGZ3QI8AOx1d13hEFkjylV0LHVz0krCeM/OAY5QWO4ZGZ9iix7cUVXVcHf3rJndAzwOpICH3P20mR0Bht39GPAZoAf4qpkB/MTd99Vx3CLSYJUqOtbj5qQ9OwcU5jWIdBOTux8HjofeOxz4+ZdiHpeINKmF2fr3fjKOAa+9rnOxouNcLr9Ye103JzWWSv6KSGTB/eZ5d/LunL88w8TM/OKumJlsrqabk6Q+VH5ARCoqt66+8DCNhfrr2/t7llR01Lp4YyncRaSsSuvqCxdNzWE2m2dqLquKjk1E4S4i14iyrr69vweAlydmMDdVdGwyCncRWSI4W8+7Y7C4Xz34pCR3X7xgWq02jKw+hbuIAKVn61pXb10KdxEpO1vf0JVmfHpe6+otSOEuskZF2QXzylyO113XpXX1FqRwF1lDFgL9uZevMjmb4/p1GTat66i4C0br6q1J4S6yRgSXXmbm8+TduTg5T0c6pV0wCaRwF0mwcksvc7k8KTOcwoVS7YJJHoW7SMJEWXpZWFe3NpjL5VnfldEumIRRuIskQKlAr7T0sjBTJw+ZNtMumARSuIu0oOByS097iouvzLG+a2mg59zJtJVeeuntTLMp28741DzdHWmtqyeQwl2kRZRbbjkzNkk256zrSC9ZS/e8425ll15u7OvhPgV6YincRZpUlNl5RzpFLu+0Fe8gDa6lA+RxLb2sUQp3kSZS6+x8IdDnc3nmcnled13X4lp6Z7qN9V0ZLb2sUQp3kVUWnJEPbuzmbduv59tnL1W8GFpudr4Q6C9enibdZtespWvpZe1SuIvUQZQA37Sug+cvTPLdFy4x0Nte8WJoudl5ps1Ip4wN3Rn6ezq0li6LFO4iMSm3pFIuwDvSKSZmsrQZXJ3OVrwYWml2PtDbyX/+1V0Kc1lC4S5SRfjCppkxMZuNvKRSLsDHJmaZy+Vps0KAV7oYqtm51ErhLmtSuWWTcIAHd6mkDM6MvQLA5g2dkZdUygX4wuuF34M3FpW6GKrZudRC4S4tL8rMOnhsdGKm7LJJe6ptSYC/cGlqcZfKhck5UmZgcGFyDiDSkkq5AM8Ul1jGJudY35XWxVCJlcJdGqrWYF6YTdc6sw6GdqrNyi6bAEsCPLhLZSHAKc7CgUhLKtUC/P23Fv6eI+NTCnSJjcI9YcqFZdTgjBqqcfSxnGDe2J3mpSszi+2izqzh1dCey+bpSLeVXDZZ0q4Y2gu7VBYfYgG0pwopHmVJJUqAf7z+Xw1ZY8zdqzcy2wt8DkgBD7r7faHjHcCXgLcAF4HfcPcXKvU5NDTkw8PDNQ94OTO9ZgmzevcRDMtsLs+LlwshuLE7zfhUIeA2b+hkNptnbHJuMTgX2lU6Vo8+Xp6YJZtztmzsYqz4MwbpNgNeDU5g8dh8Lk+mrW2x3VwujwGZ4tLHQjDn8oXvdZtB8cfFY7PzhX6trdAuuGwS/Kx0m9HX07G4S+U16zvK/j03revgwuQs41Pz9Ham2THQqxuGpC7M7JS7D1VtVy3czSwFPAf8MjACPA0ccPdnAm3+NXCzu3/EzO4E/pm7/0alfpcT7sGHDQTDa7lhs9phVu8+yoVlOBChdHDWEqpx9LGcYJ6Zz9OZbltstzCzdljyeLhK48CdfPE7lWkz1ndlyp7TdKqNK9Pz9Pd0MDmbZV3xH9PJ2SxbAv/QjoxPqUyurIo4w/1twO+5+7uKrz8B4O7/JdDm8WKbb5tZGngJ6PcKnS8n3A8cPcnoxAzd7WnOFm/HXknYrHaY1buPcmEZDkQoHZy1hGocfSwnmMPnI+rMutQ/psFZdjCkwwGuwJZmEjXco6y5bwbOBV6PAG8t18bds2Z2BdgEXAgN6iBwEGDr1q0RPnqpc+NTbOjKAFS8uAWs+FjewYyW6qPc+nCbwcI/s+G1YqDsOnLwWD36CAZzX097IXwdXru+YzGY13elXw1mh03rMoX/yRTbBfd/T85meX3/usVgDq91h49p3VuSLEq4W4n3wjPyKG1w96PAUSjM3CN89hKDG7sXZ+6VLm5BtLBZ7TCrdx/lwjIciOWCs5ZQjaOP5Qbz6/s7Fn+Osv9boS1rUZRwHwEGA6+3AOfLtBkpLstcB1yKZYQBh3Zv5/Cx00zNZSPP9JopzOrdR6WwDAZitRltlFCNow8Fs0j9RFlzT1O4oPpO4EUKF1Tf7+6nA20+CvyTwAXVO9z91yv1u9LdMqXWRiutmy7nWCv2ofVhkWSL7YJqsbN3A39IYSvkQ+7++2Z2BBh292Nm1gl8GbiFwoz9Tnc/W6nP5Ya7iMhaFucFVdz9OHA89N7hwM8zwPtqHaSIiNRHW6MHICIi8VO4i4gkkMJdRCSBFO4iIgkUabdMXT7YbAz4h4Z8eLz6CN2Ju8bpfCyl87GUzse1aj0nP+Pu/dUaNSzck8LMhqNsS1ordD6W0vlYSufjWvU6J1qWERFJIIW7iEgCKdxX7mijB9BkdD6W0vlYSufjWnU5J1pzFxFJIM3cRUQSSOEuIpJACvcamFnKzL5vZl8vcazDzP7MzM6Y2XfMbNvqj3B1VTkfHzSzMTP7QfHXhxsxxtViZi+Y2f8r/l2vKXdqBf+t+P34oZm9uRHjXC0RzsceM7sS+H4cLtVPUpjZBjN7zMyeNbMfFR9fGjwe+/cjUlVIWfRvgB8B60scuxsYd/fXF2vafxqo+JDwBKh0PgD+zN3vWcXxNNovuHu5m1FuB3YUf70V+ALXPq4yaSqdD4Cn3P09qzaaxvoc8E13f6+ZtQPdoeOxfz80c4/IzLYAvwo8WKbJfuBPij8/BrzTzEo9fjARIpwPWWo/8CUvOAlsMLMbGj0oqT8zWw/sBv4HgLvPufvlULPYvx8K9+j+EPgPQL7M8SUPCQcWHhKeVNXOB8A/L/4X8zEzG6zQLgkc+EszO1V8EHxYqQfNb16VkTVGtfMB8DYz+1sz+4aZ3bSag1tl24Ex4I+Ly5gPmtm6UJvYvx8K9wjM7D3AqLufqtSsxHuJ3Gca8Xz8H2Cbu98M/BWv/q8mqd7u7m+m8N/rj5rZ7tDxNfP9KKp2Pr5HoUbKG4E/Av5itQe4itLAm4EvuPstwCvAvaE2sX8/FO7RvB3YZ2YvAI8Cv2hm/zPUZvFB4vV8SHiTqHo+3P2iu88WX/534C2rO8TV5e7ni7+PAn8O3BpqEuVB84lR7Xy4+1V3nyz+fBzImFnfqg90dYwAI+7+neLrxyiEfbhNrN8PhXsE7v4Jd9/i7tuAO4En3P2uULNjwL8s/vzeYptEzsyinI/QeuE+ChdeE8nM1plZ78LPwK8Afxdqdgz4zeKuiNuAK+7+01Ue6qqIcj7M7LUL16TM7FYKWXRxtce6Gtz9JeCcmf2j4lvvBJ4JNYv9+6HdMisQfEg4hYslXzazMxQfEt7QwTVA6Hx83Mz2AVkK5+ODjRxbnb0G+PNiVqWBP3X3b5rZRwDc/YsUnkH8buAMMAV8qEFjXQ1Rzsd7gX9lZllgGrgzqZOhoo8B/6u4U+Ys8KF6fz9UfkBEJIG0LCMikkAKdxGRBFK4i4gkkMJdRCSBFO4iIjUws4fMbNTMwttdl9vfN83scrgAn5ndWCxC+PfFooTttfSrcBcRqc3DwN4Y+/sM8IES738a+Ky77wDGKRQnjEzhLiJSA3c/QejuczP72eIM/JSZPWVmO2vo71vARKg/A36Rwt2sUCjf8Wu1jFM3MYmIrNxR4CPu/vdm9lbg8xTCebk2AZeLRQhhGYXEFO4iIitgZj3AzwNfDVT57igeuwM4UuKPveju76rUbYn3arrjVOEuIrIybRRm2W8KH3D3rwFfW0afFyjUdE8XZ+81FxLTmruIyAq4+1XgeTN7Hyw+Mu+NK+zTgb+mUIMHCkUJ/3ctfai2jIhIDczsEWAP0Ae8DHwSeILCo/FuADLAo+5eajmmVH9PATuBHgqVMe9298fNbDuFktrXA98H7gqU0a7er8JdRCR5tCwjIpJACncRkQRSuIuIJJDCXUQkgRTuIiIJpHAXEUkghbuISAL9fwx0BJpSzVc2AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "R=[i*Ao/40 for i in range(150,240)]\n",
    "Y=[Coul_Buck(r, (r,0,0), Cl,Na,NaCl)[0] + \n",
    "   Coul_Buck(Na.latt_const-r, (Na.latt_const-r,0,0), Cl,Na,NaCl)[0] +\n",
    "   Coul_Buck(Na.latt_const*1.5-r, (Na.latt_const*1.5-r,0,0), Na,Na,NaCl)[0] for r in R]\n",
    "sns.regplot(x=R,y=Y,fit_reg=False)"
   ]
  },
  {
   "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.6.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}