{
"cells": [
{
"cell_type": "markdown",
"source": [
"# Anyonic models"
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"We solve the almost-bosonic anyon model of https://arxiv.org/pdf/1901.10739.pdf"
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Iter Function value Gradient norm \n",
" 0 8.213711e+01 1.445182e+01\n",
" * time: 0.003056049346923828\n",
" 1 6.476434e+01 1.063904e+01\n",
" * time: 0.009303092956542969\n",
" 2 5.922706e+01 1.520632e+01\n",
" * time: 0.02283191680908203\n",
" 3 4.232887e+01 1.008649e+01\n",
" * time: 0.04101705551147461\n",
" 4 3.165956e+01 8.433190e+00\n",
" * time: 0.059565067291259766\n",
" 5 2.372372e+01 5.983930e+00\n",
" * time: 0.07554292678833008\n",
" 6 2.079007e+01 5.973274e+00\n",
" * time: 0.15797710418701172\n",
" 7 1.064214e+01 2.531650e+00\n",
" * time: 0.17158198356628418\n",
" 8 9.439508e+00 2.659705e+00\n",
" * time: 0.184906005859375\n",
" 9 8.799953e+00 2.233042e+00\n",
" * time: 0.19839096069335938\n",
" 10 8.616953e+00 4.066396e+00\n",
" * time: 0.2091219425201416\n",
" 11 8.320228e+00 1.858970e+00\n",
" * time: 0.22004199028015137\n",
" 12 7.994988e+00 2.428314e+00\n",
" * time: 0.2304069995880127\n",
" 13 7.874427e+00 3.489153e+00\n",
" * time: 0.24109506607055664\n",
" 14 7.480597e+00 2.852043e+00\n",
" * time: 0.2545750141143799\n",
" 15 7.248406e+00 2.867995e+00\n",
" * time: 0.26534295082092285\n",
" 16 6.846120e+00 1.830237e+00\n",
" * time: 0.27594590187072754\n",
" 17 6.290937e+00 2.342350e+00\n",
" * time: 0.28664612770080566\n",
" 18 5.947970e+00 9.987168e-01\n",
" * time: 0.2973921298980713\n",
" 19 5.808335e+00 9.009823e-01\n",
" * time: 0.30821704864501953\n",
" 20 5.754029e+00 1.155335e+00\n",
" * time: 0.3584280014038086\n",
" 21 5.707907e+00 8.853206e-01\n",
" * time: 0.36942601203918457\n",
" 22 5.676146e+00 7.153643e-01\n",
" * time: 0.379986047744751\n",
" 23 5.660713e+00 7.717801e-01\n",
" * time: 0.388106107711792\n",
" 24 5.639805e+00 6.396866e-01\n",
" * time: 0.39899396896362305\n",
" 25 5.626301e+00 4.477176e-01\n",
" * time: 0.4097900390625\n",
" 26 5.606456e+00 3.757860e-01\n",
" * time: 0.4179420471191406\n",
" 27 5.600886e+00 6.250544e-01\n",
" * time: 0.4260580539703369\n",
" 28 5.594316e+00 3.352632e-01\n",
" * time: 0.4340829849243164\n",
" 29 5.586218e+00 3.295751e-01\n",
" * time: 0.44203996658325195\n",
" 30 5.580040e+00 3.081015e-01\n",
" * time: 0.45023608207702637\n",
" 31 5.573570e+00 1.853913e-01\n",
" * time: 0.45848989486694336\n",
" 32 5.571455e+00 1.758360e-01\n",
" * time: 0.4695320129394531\n",
" 33 5.571234e+00 3.163240e-01\n",
" * time: 0.4776020050048828\n",
" 34 5.571085e+00 2.307208e-01\n",
" * time: 0.48577189445495605\n",
" 35 5.568023e+00 2.189568e-01\n",
" * time: 0.493786096572876\n",
" 36 5.566746e+00 1.704337e-01\n",
" * time: 0.5016629695892334\n",
" 37 5.565481e+00 1.742482e-01\n",
" * time: 0.5097188949584961\n",
" 38 5.564505e+00 1.223883e-01\n",
" * time: 0.5484020709991455\n",
" 39 5.563004e+00 7.657680e-02\n",
" * time: 0.5595691204071045\n",
" 40 5.562733e+00 2.046510e-01\n",
" * time: 0.5678369998931885\n",
" 41 5.562562e+00 1.084981e-01\n",
" * time: 0.5762720108032227\n",
" 42 5.562353e+00 9.420440e-02\n",
" * time: 0.5843679904937744\n",
" 43 5.562093e+00 1.119612e-01\n",
" * time: 0.592566967010498\n",
" 44 5.561754e+00 5.788977e-02\n",
" * time: 0.6031959056854248\n",
" 45 5.561358e+00 6.454494e-02\n",
" * time: 0.6113770008087158\n",
" 46 5.561122e+00 7.479902e-02\n",
" * time: 0.619175910949707\n",
" 47 5.561024e+00 5.722902e-02\n",
" * time: 0.6274199485778809\n",
" 48 5.560832e+00 3.946113e-02\n",
" * time: 0.6353919506072998\n",
" 49 5.560752e+00 6.592560e-02\n",
" * time: 0.6435000896453857\n",
" 50 5.560708e+00 4.341524e-02\n",
" * time: 0.6515789031982422\n",
" 51 5.560663e+00 3.429373e-02\n",
" * time: 0.6621990203857422\n",
" 52 5.560629e+00 3.657749e-02\n",
" * time: 0.6702830791473389\n",
" 53 5.560584e+00 2.766974e-02\n",
" * time: 0.6783099174499512\n",
" 54 5.560578e+00 3.066428e-02\n",
" * time: 0.686068058013916\n",
" 55 5.560543e+00 2.145540e-02\n",
" * time: 0.6938381195068359\n",
" 56 5.560520e+00 2.324216e-02\n",
" * time: 0.7186830043792725\n",
" 57 5.560501e+00 1.548501e-02\n",
" * time: 0.7294111251831055\n",
" 58 5.560488e+00 1.023482e-02\n",
" * time: 0.7400960922241211\n",
" 59 5.560476e+00 7.974800e-03\n",
" * time: 0.7480790615081787\n",
" 60 5.560471e+00 9.991263e-03\n",
" * time: 0.7561631202697754\n",
" 61 5.560468e+00 4.694532e-03\n",
" * time: 0.7670650482177734\n",
" 62 5.560466e+00 4.731802e-03\n",
" * time: 0.7751379013061523\n",
" 63 5.560465e+00 2.575454e-03\n",
" * time: 0.7858200073242188\n",
" 64 5.560464e+00 3.776137e-03\n",
" * time: 0.7963221073150635\n",
" 65 5.560464e+00 2.360424e-03\n",
" * time: 0.8070690631866455\n",
" 66 5.560463e+00 2.003373e-03\n",
" * time: 0.8152060508728027\n",
" 67 5.560463e+00 1.687464e-03\n",
" * time: 0.8233280181884766\n",
" 68 5.560463e+00 1.488161e-03\n",
" * time: 0.8314290046691895\n",
" 69 5.560463e+00 1.245889e-03\n",
" * time: 0.8394739627838135\n",
" 70 5.560463e+00 1.098997e-03\n",
" * time: 0.8500490188598633\n",
" 71 5.560463e+00 1.246140e-03\n",
" * time: 0.8581669330596924\n",
" 72 5.560463e+00 1.354342e-03\n",
" * time: 0.8662550449371338\n",
" 73 5.560463e+00 8.755295e-04\n",
" * time: 0.8743751049041748\n",
" 74 5.560463e+00 8.394436e-04\n",
" * time: 0.8961079120635986\n",
" 75 5.560463e+00 4.868165e-04\n",
" * time: 0.9068179130554199\n",
" 76 5.560463e+00 3.382988e-04\n",
" * time: 0.9175510406494141\n",
" 77 5.560463e+00 3.029369e-04\n",
" * time: 0.9284369945526123\n",
" 78 5.560463e+00 5.095378e-04\n",
" * time: 0.936439037322998\n",
" 79 5.560463e+00 2.615929e-04\n",
" * time: 0.9474270343780518\n",
" 80 5.560463e+00 2.313423e-04\n",
" * time: 0.9580779075622559\n",
" 81 5.560463e+00 2.405130e-04\n",
" * time: 0.968635082244873\n",
" 82 5.560463e+00 1.953079e-04\n",
" * time: 0.9794259071350098\n",
" 83 5.560463e+00 1.951937e-04\n",
" * time: 0.992811918258667\n",
" 84 5.560463e+00 1.659792e-04\n",
" * time: 1.006005048751831\n",
" 85 5.560463e+00 1.659795e-04\n",
" * time: 1.0391590595245361\n",
"e(1,1) / (2π)= 1.7391793785607679\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3deXwURd748eo5cnEkhCREAhIhnFECixy7oEYgUTkkCqJZXQHRR3B3fVj34WF1FVnvrD7rhayuIiyIyKKiLwRdBQUE5RJFDkHkJiQcIjdJZjLz+6Pd+U1Xh5rM0KEzzOf9yh/T3TXdlc4klapvf6s0v98vAACIVQ67KwAAgJ1oCAEAMY2GEAAQ02gIAQAxjYYQABDTaAgBADGNhhAAENNoCAEAMY2GEAAQ02gIAQAxzXV+LuPz+R6c9MiEP//5/Fwu2nm9XpfrPP1oLhjV1dVOp9PuWkST6upqh8OhaZrdFYkmMf67mRzvDllm69YjZ854wzrtrl0/HDiw4e677460XufqPP1EKysryz3a4l2Hz8/lAACWu7H9RSHL3Hzz/PXrD4V54o1FRcdsbAgZGgUAxLTY7eMDACynaSLqhtvpEQIAYho9QgCAZTRNC/cJLL/f5i4kDSEAwFIMjQIAEEXoEQIArBOFD8vQEAIALBNBjFDYHSNkaBQAENNoCAEAMY2hUQCAZSJIqLd7ZJSGEABgnQhihLY/XMPQKAAgptEjBABYR4u+hHoaQgCAZbTwhzptHxqlIQQAWCaSPEK7W0JihACAmEaPEABgnQhihKHK//jjj2+88cbRo0evv/76rl27mgv4fL5//etfmzZtysvLGzp0qN4lff/99w8cOBAok5GRUVRUVOP56RECAKykhf+lcOLEiR49eqxevdrhcPTr1++TTz4xlxkzZszTTz/dqFGjRx999A9/+IO+c+vWrV/9xxNPPPHuu++e7RL0CAEA9desWbMuuuiiWbNmCSFSU1OfeOKJgoKC4AJ79+6dMWPGrl27MjMzb7nllvbt2z/wwAMZGRn/+7//qxeoqKho3rz5HXfccbZL0CMEAFhG035+XiYMyj7hp59+et111+mvr7vuumXLlnk8nuACS5Ys6dy5c2ZmphDi4osvzsnJWb58eXCBd999NyUl5aqrrjrbJWgIAQDWiWBgVDk2un///oyMDP11s2bNfD5feXl5cIHy8vJmzZoFNjMzM/fv3x9cYOrUqaNHj1Y8y8rQKADAMnqHMLy3aNqaNWuGDx8evPPaa6/VBzOdTqfP59N36i+cTmdwSYfDESgghKiurg4usGvXrmXLlk2fPl1RARpCAIDNmjdvftNNNwXv6dChQ+BQWVmZ/nr//v0ulyvQQQwUCO4ClpWVNW/ePLD5+uuvFxYWtmzZUnF1GkIAgGUiWH1CCJGVlSU1hAEDBgx49tlnH3roIafTOW/evMLCQpfLJYT47rvvUlNTmzVr1q9fvzvuuGPHjh2tW7f+7rvv9u7dm5+fr7/X5/P985///Nvf/qa+OjFCAIB1rI4R3nTTTQ6Ho7Cw8O677/7rX/86ceJEff/IkSNnz54thMjIyLjvvvsKCgrGjRt33XXX3X///cnJyXqZjz/++NSpU4MGDVJXmR4hAKD+SkhIWL58+YIFC44fP/7ggw8GBjlfeumlwDMyjz/++IABAzZt2lRcXNyzZ8/Ae3NychYvXhwfH6++BA0hAMAymgiRDlHTW0JISEgYOnSotPPyyy8P3uzdu3fv3r2lMjk5ObWpAA0hAMA64ccI7Z5zm4YQAGChKFyPkIdlAAAxjR4hAMAyESTU2z42SkMIALBMBHmEdreDDI0CAGIbPUIAgHWi8GEZGkIAgGV+XlgpzLfUUWVqiYYQAGCZKOwQEiMEAMQ2eoQAAOtE0CW0uwtJQwgAsE4EMcI6qkmtMTQKAIhp9AgBAJbRwk+Qt/uhURpCAIB1NC38dAi7W0IaQgCAdaLwYRlihACAmEaPEABgmSjsENIQAgCsE8kUa3Y3hQyNAgBiGj1CAIB1onBslIYQAGCZCPIIbW8IGRoFAMQ0eoQAAOtE4VyjNIQAAOsQIwQAxDJNi765RokRAgBiGj1CAIBlNKGFmyAfNQn1p06devDBB/v06ZOfn//0009XV1fr+0tLS2+//fbLL7989OjRBw8erLN6AgCigRb+l91q2xD+6U9/WrZs2UsvvVRSUvL6668///zz+v4bbrihcePG06dPF0IUFxfXUS0BAKgjtR0aXbt27V133ZWXlyeEGD58+Nq1a4UQq1at2rZt24oVK9xu9wsvvJCenr558+ZOnTrVYX0BAPVYBA/L2N4prG2PsKioaPbs2bt27fruu+/ef//9oqIiIcS3337btWtXt9sthGjQoEFubu769evrsLIAgPpN+3na7fDYW+faNoR33313RUVFt27devbsmZOTM2TIECHEwYMHU1JSAmWaNGlytjCh3+8PhBUBANHI4/GELnQBxwiHDRvWq1evQ4cOHT582OVy/e53vxNCNG7c+PTp04EyJ0+eTE5OrvHtmqY5nc5zry4AwC76+N+Fp1YNod/vX7p06a9//WuHwxEXF3fLLbcsWbJECJGdnb1t2za9jM/n27FjR3Z2dp1VFQBQ70UyNmpzlWvVEGqa1r59+4ULFwohfD7fwoULO3bsKIQoLCw8ceLEggULhBBz586Nj4+/4oor6rS6AID6LApHRmv91Ohrr732m9/8Zvr06ZWVlRkZGW+99ZYQIj4+ftq0aaNGjUpJSTl58uQbb7zB+CcAILrUtiHs1avXtm3bDh065Ha7gx+QGThwYGlpaXl5+UUXXeRyMU8NAMS2Oph0+8CBA9OmTTtx4sTgwYN79eplLuD1emfOnLl169bc3Nxbb73V4XAE9s+ZM2fDhg1NmjQZOnRoTk5OjecPb67R9PT04FZQ53a7W7ZsSSsIAIgoeULVEh4/frxHjx7btm1LTk4eOHCgHqSTjB49+pVXXsnKynr++ef1ZzmFEFVVVQUFBVOmTElJSfnpp5+++OKLs12C1gsAYJkIVqhXl58xY0Z2dvbUqVOFEMnJyU8++eSAAQOCC+zatWvOnDl79+5NT0+/8cYbc3JyJk6cmJmZOWXKlJMnT3755Zch+2msPgEAqL+WLl16zTXX6K8LCwtXrFhRVVUVXODzzz/Py8tLT08XQmRlZbVr127FihVCiAULFowYMeLDDz984YUXvv76a8UlaAgBANaJ4LFRZY+wrKwsIyNDf92sWTO/319eXn62AnqZ/fv3CyF27tw5efLkd9999+DBg9dcc40+J3aNGBoFAFgmksRATVuzZs3w4cOD9xUUFNx1111CCJfLFZiYzOv1ClNev8vl8vl8gU2Px6MX0DQtLy9v2rRpQoi8vLw//vGPI0eOrPH6NIQAAJs1b978pptuCt6Tm5sbOKT38IQQpaWlLpcruP+nFygtLQ1s7t+/v3nz5kKIrKyswEkuvfTS0tJSr9dbY7yQhhAAYJkIVp/QhMjKypIawoBBgwaVlJQ89NBDLpfrnXfeGTBggJ6w/s0336Snp2dlZfXv3/+OO+74/vvv27Vrt2HDhrKysvz8fCHEkCFDPv30U/0kq1evzsnJOdtTMzSEAID6a+jQoZMnT87Pz2/Tps2CBQv+/e9/6/vvvvvu4uLicePGpaWl/elPfyooKLj22msXLlw4ceLExo0bCyFGjx49ffr06667rmXLlvPmzZsxY8bZLkFDCACwTvgxQnXx+Pj4JUuWLF68+KeffiopKcnMzNT3v/rqq/qTokKIiRMnDhw4cPPmzWPHju3SpYu+s3HjxqtXr160aNGpU6cmTZqkj5fWiIYQAFCvxcXFXXfdddLOzp07B29269atW7duUpn4+PiBAweGPD8NIQDAMtG4Qj0NIQDAUnY3bOEioR4AENPoEQIALBNBQr3tK/PSEAIALBNJHqHdQ6k0hAAAy2gixLJKNb/JVsQIAQAxjR4hAMA6EaxQbzcaQgCAdYgRAvWdX7FVr8l/K+z+2wFcMGgIAQCW0cJPh6BHCAC4gBAjBGwQ3min8vg5vTk8of9WGEvIl/arTsA4KuwSjXmEpE8AAGIaPUIAgGUimGLN9i4hPUIAQEyjR4go4A8RmlMFCcN6a+jjoapSe37Tf8Hytj+Mw9JBOYCoDCgK+/8jB+xEQwgAsEwkD8vUTU1qj4YQAGAZTUTfMkzECAEAMY0eIeqB0FG9cKKA0lHlYfmo6dr+cEKOIYJ6Ek0+l7x4jbSlDikqL6ZJVTMV9huDiKQhInIRJNTb/QGjIQQAWCf8GKHtaAgBAJaJII+QGCEAAHaiRwgbhEjGU4f1TNt+vyrsF95Rn7JipvLhzT0qx/zM/wUbUwOl8g5jGE+TCktHpUuHildKoVOpvBRhjLaBL0CNhhAAYJmIhkbrqC61xdAoACCm0SNE3QgxghhmDoN6eNNn2Pb5VIXlo8pNIYS8Q1kTtZDDlQ71WKhx06EcKZWPGv/jreEfcHmXapCWXAsoRDCzjO2fGRpCAICl7G7YwkVDCACwjKaFnQ5B+gQAAHaiRwiLhJj2zK84qo75mff41JvV4RxVFjbvkWOK5xAjlGJ+whTYkzedyk11YeWZhTmIaCpgOKrMtZBnjou2UTKcqwhWn7D7Q0JDCACwTgRzjdqNoVEAQEyjRwgAsIwmNHkdldq8yVY0hIhUWEHBcKY9k+J2whSoq642zIRW7TUe9YZz1Hgq86Wrq5UxQp8q9innDSrzAoUQTimS5zJsOp2G8RunS9rUwjjqlIeCHH7VpjxyJAUFwwoZCtv/6KFuRbJCvd0fCRpCAICloqxDSEMIAKjf9u3b99prrx09erSoqCg/P99coKqqaurUqZs3b+7cufOoUaNcLpcQYtWqVUuWLAmUGTNmTHJyco3n52EZAIBltAgoT3j06NGePXsePnw4Jyfnpptuev/9981lRowY8dZbb3Xu3HnatGljxozRdy5ZsmT27Nk//YfPd9bFZegRonbUs4WK8IKCYeX2CVNgz2sM+3k9xqOeasNmleGoRy4sBRTlXxV1EFFatskvTdEppKCgobCU6idCRQFdbsOm27jpilMVdrmdhnq65DssXcsvpSEKKQ1RencYIUNBouGFzvK5Rv/5z3926NBh8uTJQojExMSnnnpqyJAhwQW2b9/+3nvvlZaWpqamDh48ODs7+5FHHmnevLkQokePHk899VTI69MjBADUX59//nn//v311/369Vu1alVVVVVwgRUrVnTp0iU1NVUIkZmZ2aFDhy+++EI/tGXLlscee2zatGlHjx5VXIKGEABgHS2ir7MrKytLT0/XX2dkZPj9/rKysrMV0MvoBZo1a9apU6fq6uo5c+Z07Nhx9+7dZ7sEQ6OolRqmErNuLFTKcPCaxielAUxPpWHw01Nl2KySjlZKR41Do8b3SgOnIlQyhil9Qj00KmU4mIZGXerBT8PwZly8sXC8dNSw6Y4z3n/jUSGENFbq8ofxL3JYI6XCnF9R+yshKkSwMK/QVqxYUVBQELzz+uuv//3vfy+EcLvdXq9X3+nxeIQQ8fHxwSXj4uICBfQycXFxQoiRI0eOHDlS3zls2LCSkpIpU6bUWAEaQgCAZSKYYU0Tom3bthMmTAje2bp1a/1FVlZWaWmp/rq0tNTtdgf3/4QQzZs337dvX2Bz3759WVlZ0iV69uz52Wefna0CDI0CAGyWkZHR3yjQEA4ZMmTevHl6XHDOnDmDBg1yOp1CiNWrV+/Zs0cIUVhYuH379s2bNwsh1q1bd/DgQT3F4tSpU/oZPB7PggULLr300rNdnR4hAMAytUiIqOE9ioM33HDD3//+9z59+mRnZy9ZsmTRokX6/t///vfFxcXjxo1r0qTJxIkTCwoK+vfv//HHHz/22GMNGzYUQlx55ZVNmjRJS0tbu3ZtkyZN7r///rNdgoYQNTMtMBRi+5yCgh4p5icH6qSwX1WFYbOywmvYPCNtSu81HJUDiqZLq/Mr5GWblIkBDjlGKI/HyAkSxowIKewXl2D45Y1PNBz1ykdDrDbl9xmjhvKPX44pKsghQ3P6hJRf4TeGDIkZXgAsTZ9wu92LFi1aunTp8ePHX375Zf3pUCHEjBkzAq8nTJgwaNCgzZs333///R06dNB3fvzxx1999dWJEyfuu+++7t27K5pnGkIAQL3mdDr79u0r7Wzfvn3wZm5ubm5ubvCepk2bFhYW1ub8NIQAAMtEMul23dSk9mgIAQCW0TQRdvqE3S0hDSH+o4ZUwaCD5inWpOiYMbgWVlBQyu2TwnjCFParkDZPGTbPnDYePe0xnkoVX5SSDoUpalinMUIpKCilBsbLYT9DzRMq3cGb3iT1alPyL76cEOlX/2VQhQzltehC/Y2T/wgyAVu0iyx/wlakTwAAYho9QgCAZaKwQ0hDCACwkNV5hOcBDSF+pl5WqYZlmHxhJA5K04dKgTcpKCiFAIUQZ04Z4nxnThoKnD5ZZTxqLHxKChmqkg6raogRGvZIwU6fMfYm3TPpt9shr7Ik//K73VKmoDFGmGj4bU1IMmxK0U2vxxAyrA610JUyQCw0zSVtG7ccxk3DyRymP3GavEu1bJPtfQXEgjBihN98801RUVGHDh2uvvrqFStW6Dt37do1dOjQjh073nLLLfv376+bSgIAood1S0+cH7VtCLds2VJQUJCfn//ee+9NnDgxKSlJ319UVNSmTZuFCxempaXdfPPNdVZPAEAU0PMIw/2yV22HRh955JFbb7113LhxQojABDZffPHF3r17n3jiCZfL9fTTT6elpW3YsOGyyy6rq8rCUupJ1NQzqIlQQ6PymvLGicqkEUhpfFIaCBVCnD5h2HPquGEs9JTx6OkThqOnjSOlUq6FNAwrTd4mTJkepqFRKevAuAyT8fdbWpLePDTqUg6NJhiHRisrjEOjVcb0iXDqaSbVTBrMlDalkU65sOlvnF8+u3KBeyZgizYRzDUadkzRarXtEa5Zs6Z9+/ajRo267rrrnn/++erqaiHExo0bu3Tp4nK5hBCJiYm5ubkbNmyow8oCAGC12vYI9+7d+9e//nXy5MnJycljxoz56aefJk2adOjQoeTk5ECZJk2aHDx4sMa3+/1+ve0EAEQpj8fjdrtDFKofYb+w1LZHmJycfOeddw4cOLBPnz5/+ctfZs2aJYRISUkJLPgkhDhx4kSTJk1qfLumafoKUgCAKBW6FYwoRmh7w1nbHmGbNm0Cnb+UlJSTJ08KIS655JKtW7fqO30+3w8//BBYShHRR44KGrdMQSU5KGhMJJBiVB7lOkpSoE7KjhCmoODJY6rNU1KM0BhBVM+4Zo4RepQxQtPMZKoYoSZPsRZe+oRpASljUFBeLkpVTzN1OFNK/JDmipM2NaehJlJyhTAHEUOkbiiPAlaobY9wxIgRb7755unTp6urq19//XV9RYz+/ftXVla+8847QoiZM2c2bNiwd+/edVhZAED9FoUPjda6Rzh69Og1a9a0atXK7Xbn5eVNnz5dCBEXF/fGG2/cfvvt9957r9vtfvPNNx3yopwAgFgShXOs1bYhdLlcr7322osvvqhpWkJCQmB/v3799u7de+zYsZSUFNsfgQUAIFzhTbGWmJho3ulwOM72jAzql3NIHDQHmUyJg6qFlqSJyuRllU6ppkwTpkxBU4yw0lD4uCqPUJpxzRR4M+cRSgsYqbInQ0yx5lDF4YQQLrdhQCWu0hgjrDT8tpoqZpzsTRkUNE1yZspxlDdVk8NJm3J80SHXxCHdtHDSCuXS/O9dH4WfR2j3D5K5RgEAlokk6Gf3PzQ0hAAA60RhjJBnWwAAMY0eYQwJsfSOMnHQZ44RVksxwnAmFzVG5s6clvIIzXONqjIF5alHjxunHj2pyiOUEgel9aGEaQEpUyjOUFi6R6YZOw2bUnKeEMLlMuzxVqluqTx9qFxx6dKhpjk17pFmPXXFOVSbbimCaKiK0xQK9RmjhmGlFRIirP8imWvU7p8kDSEAwDJRODLK0CgAILbRI4xdpnwJaVOVGCBCpk8Yh/U88rpLxinWTqsyHIRpeFOaNU0+KudLSHOqqfIlvKahUWnqONN6RobCfuNNlQZ8TAvWm65lHCytrjZsqlM1JOphWGkwUwjhjjOMhbrjDTctLt4491u88aYZ3+tyGzNM3KYh+RCfNNXop+1dB4QWhV1CGkIAgGW0CNYXpCEEAFwwIsgjtLsdJEYIAIht9AgvaCESJlRHzzF9wiOnTxg3jekTppChHCOU5mBTxxRNCy15FZtSvoT0XQghquUsBXWgTpVAIf2brFWb8gqcyvnb1LOmyQFIKR3COHlbvLw4aHyC4bbEJxr+MlQlSoFVw9G4BMNNc8sJJ/J/29L35ZDnWFMnoTDjWr1HjBAAENMiyCO0e8EGhkYBADGNHiEAwDJMuo36RR0SVM+pJgfDQsYIq1XzgXmNyzCpZ1yTwnhCiAo5zmc8m7SUUoXqWh7lvGXmGKE6e8+vzuaT0gqlBYZMt9TvMxQIlSlo+L6kNZ5cbumOGYKCFQnyHU6Q7rC8OpVxBSjjT1O+h8rMSxHqo+VXz7hm/IupDicCtURDCACwjBaF6xESIwQAxDR6hAAAy0SSUG/3oDYNYSwJI54VImQoTMEzn3GuUXVaoUc5E6kU5DPvCS8KqMwUrJaT3kzfpnrO1RAxQmPMT5NChqbiposrynurDNsep5S4qbxjoe6wR3lLPcrAqvRJMGeghsjFDOdTaveIGmoQyTJMdreENIQAgHpt586dL7/88rFjx4qKiq699lpzgYqKipdeemnr1q25ubljx46Ni4sLPjpnzpzDhw//9re/Pdv5iRECAKyjhf+ldOTIkV/+8pc+n6979+4jRox45513zGVuu+22jz766Kqrrpo3b96dd94ZfOjLL7/83e9+98QTTyguQY8wdplGpFRjUubpveRlmORVmVSjZ171SKlHHrjzGvdIBeRUDa/qaX5pyjR5gNc8Z5qcTxHGWkhy+oR0TJ5aTF663TgyKtfcUa36ruU7LN8x+Q5LBcJKMpHvsC/E0Khyqjj5c+iXpqlTvRX1g9V5hNOnT8/Ly3v66aeFEE6ns6SkZOjQocEFvv/++wULFpSXlycnJxcWFrZs2fLxxx9v2bKlEKKysvK3v/3txIkTn3rqKcUl6BECACyj/fy4TFhfKitWrOjbt6/+um/fvmvXrq2srAwu8OWXX3bt2jU5OVkIkZ6e3rFjx5UrV+qHJk2aNGzYsI4dO6rrTEMIAKi/ysrK0tLS9Nfp6el+v7+8vDy4QHl5eaCAECIjI6OsrEwI8c0333z88cfjx48PeQmGRgEA1olo9YkVK1YUFBQE7xs4cOC4ceOEEPHx8VVVVfpO/UVCQkJwyfj4eI/n/685U1lZmZCQ4PV677zzzsmTJ7vd7pDXpyGMIabwl19xWJ48rIbgmWqiLDmCqEyuUG/W9HZpdjfjpZVTfKmf3a8p5hdWUFBFeq8mzx4mRw3VGSy+cL5r9Q2sqYD6B6SqScjJ+YSckaJcWUkqrBEyrO+08PMCNSHatm07YcKE4J1t27bVX7Ro0WLfvn36671798bFxaWnpweXzMrK2rt3b2Bz3759WVlZ33///bfffnvbbbcJIc6cOXPw4ME2bdosWbJEjx1KaAgBADbLyMjo379/jYeKiooefPDBhx56KCEh4c033xwyZIjD4RBCLFu2rEWLFq1bty4sLBw1atT69evz8vJWrVr1448/Xn311S6Xa8uWLfoZVqxY8T//8z+ffPJJZmZmjZegIQQAWEbTwk+QV5YfMmTIq6++2qNHj+zs7DVr1ixevFjfP378+OLi4nHjxiUnJz/++OPXXHPNlVdeuXTp0pKSkqSkJCFE69at9ZI7duxwuVyBTTMaQgCAdaxeod7lci1cuHD16tVHjhx54403GjdurO+fM2eO/qSoEOK///u/Bw4cuGXLlpKSkksuuUQ6Q+/evb/44gvVJcKsL2JGqOBZqGWbVEfl7D1l9Mu8J8Q6PtKl1VPHhQyFnkNQUM18ZlOioRRLU0cQDZvy6lGh77AxxKj8AYX1sw75yQm1WhiiTCRzjYYq4HA4evXqJe3Mzs4O3szJycnJyanx7YmJia1atVKdP1QFAAC4kNEjBABYJhon3aZHCACIafQI8bNzjNSoZy6VJzJVRolqiiqpMs+U+ZDmOVTrc0hKvciQ6paql4dS30Bz+RA/oHBqEq76/OPBhYqGEABgJWsn3T4PaAgBAJaJJEZYR1WpNRpC/OwcP4uavKkpNs2lDVumqki/V6bNMGqiHm+0m/qHoLql8jcp31LVDTSXD/EDCqcm4bL9byLOldV5hOcBD8sAAGIaPUIAgGWiMX2ChhAAYJkIZpaxHQ0hziJ03M646TAGihyqow71plMesXc4pZOr3i5fWvmNSP+K+s3fpnXLMMlnNv+xCBHYU71d+q5N90S6w/K1pXuu/gGF9bMO+cmxPT4EECMEAMQ0eoQAAMtEY4yQHiEAIKbRI4whcgKdHGUyhL/k/9HMmWfqKJRx0+mSNh3KTflaprerIl5yuMtpWCVIqzZW2yd919KV5cWPNC3ykKH5FqoLqKOwob5r9f0PeYfD+AGpI7jSphDyN2b6pCkLy+dCfWTpurznAw0hAMA6UZhQT0MIALBMBDFC2xEjBADENHqEsSvUnJyG8Jc50CNHAeWgoCrI5HIbNt1xxk23U7qWy7hHKiCdzWW8VrUxeOZzGoN8PuM3Zor5+RxSUFAKWUUeJDT/0+ww7pLvsBT2c6q+a/kOy3dMvsPuOOMdjlOdTf7hOqVN1XchavosBbN25lKcfxEk1NvegaQhBABYJwpjhAyNAgBiGj3CWKJegCicZ/eFeRYu5fP3bvVYaLxhXC4uQR64k/bIm5WGTW+VIV+iuto4NOpT5z/4hHKHaeV3ubiCeso0UUPOiXJ4OU51S+OUt7SGO2wsL/1ETGPXqorJeR2mT456hrZwlqJCfaTVEGcJ+RYm3QYAXCi0CPII66YmtUdDCACwTmQxQltXyCZGCACIafQIL2TqmGBYQcEaYoTqp/mllAbj0z8yYlgAAB2eSURBVPlSRCo+wfA5jE+UP5YJxj1VFV7DZqXhqNdjCOvJQUGf6j9P85BOdbXq7X553SZV3DVk2FVOkHCdSxRQuqWGowmmOyzdc+knIl1L+mnKyRXKyd5EDes0KUOnIXItUP9EkFCvhZeIZDkaQgCAZSLJI6ybmtQeQ6MAgJhGjxAAYJmI1iOso7rUFg3hBc208JK0rSgtrTdUw0RZ6oCWWwpoGTflCJYh5peQJH8sKyuMMcJKd/CmHBSsDmOlJM1RbThVlenbrA7j5H7jLZVnCzOe2xw8k+OsyqCgFNVLSDLck8QGqqMJDUxR2CRVTDFODhlGnlYoTJ+lkKtTqY7a/QcUFwYaQgCAZSKIEdr+Dw0xQgBATKNHCACwjCbCjxHa3SWkIYxd6tCMKbtLDrWZJsY0bMozYcrxLWNamzHm56mUJ/z0VBkjeXJQ0LDpN00XGkwO1Bm/C49TfrPXG8a1lLO3ypNqOpzyeIxLnThozBSUwn5JDZWbjVRHRaiYovTzcstphVKMMORco8ZNZR6h3WNmCN8Fv/qE1+v99a9/PWbMmMCe77///tprr73kkksGDx68e/duq6sHAEDdCq8hLCkp2bBhw6pVq/RNv99/ww039OzZc+XKlR06dBg+fHgd1BAAEDX0h2XC/bJXGEOjW7ZsmTt37n333ffCCy/oe5YvX37w4MGJEyc6nc5HH300PT39m2++6dKlS91UFecqRDKFVFg5hChCpU+45PQJ40pJxkf/pYWTvB7DQKgQwuuNC95UT3smMQ1IGodw3dJsbaZLK3Mz1Is6qedUc5ryCtQ3TcphkAYzpdHOBo0Nd6xBI8NmknFTCJFofHuinE2hnHHNuN596PSJsOZUY6Q06kQwxZrdP9jaNoQ+n++uu+566aWXgsc/N23a1KVLF6fTKYRISEjo1KmTvqdOagoAqP/qIEb4/fffv/TSS0ePHr3hhhuKiorMBU6dOvX8889v3ry5c+fO9957b0JCghBi9erV77zzTmlpaWpq6rBhw6688sqznb+2Q6N/+9vfunTp0rt37+CdP/74Y6NGjQKbKSkphw4dqvHtfr/f6/XWeAgAEBVs+TN+6NCh3r17N27cuLCw8J577pk9e7a5THFx8YoVK4qKihYtWjRy5Eh954YNG5o2bTp48OCsrKxBgwZ98MEHZ7tErXqEe/bsee655xYsWLBjx46DBw9WVlbu2LHjkksuadKkycmTJwPFjh8/npqaWuMZNE1zuXhCFQCiWG3+jFveIZw2bVr37t0fffRRIYTP53vmmWeKi4uDC2zZsmXRokUHDhxo1KhR3759mzdvvnv37latWo0ePTpQZu/evR988MGgQYNqvEStGqf9+/fHx8ffeOONQoiTJ08ePXq0oKBg48aNrVu33rJli9/v1zSturp627Ztbdq0qc0JYQ9lkFBeB0UK5JjGDkKkTxjjRu44w8njEwyb1cYUBSkEKMwTm4UIChof31euD1V5xlDPqgo5RugxBiyrvepVmYxTrMn3UHXHhBBut2plJTnnJEmVICEFBRsmSyFDOX0iqaGhgDQHm7Qqk1u9DFPo9AnlMkxSaYKE0UYz/UxDv0VZfOXKlfn5+frr/Pz8ESNGVFRU6IOfgQJdu3bVhydTU1Nzc3NXrVrVqlWrQIG9e/cuX778j3/849kuUauh0V69em3/j2effbZTp07bt29PTEzs37+/z+ebNWuWEOLVV19NTU391a9+VZsTAgBQG+Xl5U2bNtVfp6Wl+f3+8vLysxXQy5SVlemvZ8+erWnaxRdffNlll912221nu0TYU6w1bNgwKytLf+1yud56662JEyempKQ899xzb775ZvgPCwEALiBaJF9Lliy53KikpEQ/X0JCQlVVlf66srJSCJGUlBR8wcTERI/HE9isqKgIFCguLvb7/du2bfvhhx8efvjhs1U57Ljd9ddff/311wc2+/Tps2PHDqmjCgCITZHlBXbp0uWZZ54J3hPocbVo0WLPnj366z179iQkJKSlpUklAwWEEHv37m3RokVwgZycnLFjxz7//POPPPJIjVe35gEWWsELQThphcIcIzSG4vwuY76dMfNMSr/zVbsUR4Up9ibXTZmfp07Oq0gw5hGaYoReZYzQtCqTKkZoyrw05xGqYoQJ0rpLyjxCKVNQCgpKWYZCiMSGxpNLiYMJyinWXKpPQg0xwnASBxF1IlqPUEtJSenWrVuNR4cOHTp+/PiHHnooKSnpjTfeKCoqcjgcQohPPvmkVatW7dq1KywsvOOOO7766qtu3botX7782LFjekzx4MGDGRkZQgifz7d48eJ27dqdrQI8yQkAqL8GDRr0+uuv/+IXv2jRosWWLVsWL16s73/wwQeLi4vbtWvXuHHjkpKSa6+9tmfPnqtWrfq///u/xMREIcSAAQN8Pl+zZs22bt3aqFGj+fPnn+0SNIQAAOtYnT/hdDrff//9b7755tixYz169NAbOSHEe++916BBA/312LFjBw8evGXLltzc3IsuukjfuXLlyo0bN/7444/Nmzfv0KGDop9KQwgAsJAW7rJKtSlvnrMs0ODpWrRoIYUGXS5XLWc6oyGMXdK/R/5w0gqFeQ5PY3F56lHjUb/Padw0Rtpqrm9wXVSxN9NiRsZsPGNQMOFMiBihOo9QXpVJec/Us7OK0HmE0kpJqqlHE+VVmeKMR+Vf/MQG0rpL6slFpcRBVVBQqyFGGEbiIM+h4zygIQQAWEaLwn9faAgBANaJwoV5aQjxM9P8a/IOubw0k5lU2jjQ6vI7pMPGrRCfQ6kqpswN5exu8YZLS4N+lRWGoVGPaRkmj3GJKGk2ODnTQ7lEvWlSOtPQqHGpdylLQZrnTD3jmrSOkpRrIWVHCPNYqDynmnGgW7nQkrRpHho15UuoR0oRZSJLn6ijytRS2DPLAABwIaFHCACwTAQzy9g+DEBDCACwDjFCRDFlTDD0v3jKVAETp+KYpskfS/XKSnKqRpxqTrWqRGO+RKXhWlJEUAjh9Zy/GKGUluBWfiNxypChKaComjLNfHI5KGgMu8oLLUlBQWndJHOIUP3RsPtvImIQDSEAwDIRPCxj+38/NIQAAMtE4cgoT40CAGIbPULUTD0Bmwg1B5sjvH+xjCGrGqZzk06uWs/ItO6SITVQCgp6qwxHPR45RigFBaUp1kyTwxmXYZJmC1NWW5iihm63FPs0xgjjVUmHppifah0lYbppUqagOigoz6kWOi+QSdQudNH2M6UhBABYJwoT6mkIAQCWiSSP0O4eJDFCAEBMo0eIWjH/x+aXQ0GqHLqwQoaaJpeWpqsMK0bolUJrCYaYnzpN0LxHShz0G4uHihEaCjtMk3A65ciclMynjiCqCkuJgK5Q4Un19KEhgoLKqURF9MWPEKYIHhu1Gw0hAMAymhZ+zM/usVGGRgEAMY0eIWrHPFGWNLnYOYyUmsbWQiz5pDkNI5Kmld8NR11uw9ncUjqEcYl5n1e+dLU8FqpOnzBWW9oMtXS7U8oSMQ5gOpUTy0njwyGOOk2Dz8oEFS3Ez0uIs2/X8K9+tI2bISwRLMxrd4eQhhAAYB1NaFqY/+zY3Q7SEAIALBSFc6wRIwQAxDR6hIiUtNpOWCFDeQI3w5b5vzPTyj7GZ/0dhvc7jSFDnzFG6Kt2GDeNR31yjNCUL2HclIuryN+FKUaoDtSFldIQ4qjp0ur4pVxz+c3KoKDd/+zjfAs/od52NIQAAMtEsAyTZn5A7vxiaBQAENPoEQIArBOFD8vQEMIi5xIylM5kimD5ldE1hzFu53OownghYn7mGKG0w29ZjNAcSHEov031WlTSyUyJgCFWPpLHsuSgrFz87GXt/6MGe0Uy6bYpB/c8oyEEAFgmgjxC2/97IkYIAIhp9AgBAJaKsg4hDSHqiDLzzO9Xh6xM8QIpxCgHII1ljfEwvzGO5zCWDjcvUDpbeJENOfAW4rdfnXdoyq0M42gNM8eqrx3GQcQ6FuYFACDK0CMEAFgmkoT6OqpKrdEQwgbyr4ky10KESrcIMXDqV13ML6UsmK7jN1VOURH1MkyyGsYnzyGHQc7NUF4q1NpIZEQgcuQRAgBiWgRzjdrdEBIjBADENHqEAADLaLV4Itr0Fpu7hDSEqAdCPNtvivPJW+oIomo75KT3fun9YU2qpqYO+tVQPpyD4Z7b7rEpQGHjxo0vvPDCsWPHioqKiouLzQWOHz/+17/+devWrbm5uePHj2/QoIEQYtWqVfPmzdu5c2dGRsadd96Zl5d3tvMzNAoAqL8OHDhw5ZVXZmdnDx8+fMKECTNmzDCXufnmmzdv3jxq1Ki1a9eOGDFC3/noo48mJCQMHz48NTX1l7/85fr16892CXqEAADLRLYeoeLo66+/3qdPnwceeEAIUVVVVVJScvvttwcX2LRp07Jlyw4ePNigQYNf/epXmZmZO3bsaN269fz58/UzDx069Jtvvpk3b97ZOoX0CAEAltFnlgnrSz1yv3r16iuuuEJ/fcUVV3z77bcVFRVSga5du+rDoSkpKZdeeunatWtFUPvq9/v37dt30UUXne0S9AgRBUL8fynNmhYi4uhXHTZP7iZfyrJgWugTnUPMkERA2CaCPEKl8vLypk2b6q/T0tL8fn95eXl2dnagwIEDBwIFhBBNmzYtKysLPsPkyZNPnDgh9SOD0RACAGy2ZMmSyy+/PHjPjTfeqA+HJiUlBbqAZ86cEULonb+AxMTEysrKwGZFRUVwgTlz5jz55JOfffZZYmLi2a5OQwgAsExEMULRpUuXZ555JnjnxRdfrL9o2bLl7t279de7d+9OTEwM7v9JBfQyLVq00F/Pmzdv3Lhx//73v9u3b6+oADFCAIB1wo8RappISUnpZpSenq6fb9iwYW+//fbJkyeFENOnTx86dKjD4RBCfPDBB5s3bxZCFBYWlpWVrVy5Ugjx2WefnTp16uqrrxZCfPTRR2PGjJk/f37nzp3VVaZHiOgXKg3RQAooKk8V6lLnFWE/xKYBAwbMnDmzc+fOWVlZ+/btW7x4sb7/0UcfLS4u7tSpU8OGDZ999tlBgwZ17dr166+/fvHFF+Pj44UQ99xzT2Vl5c0336yXHz58+JNPPlnjJWgIAQD1l8PhmDNnztatW48cOdKtW7e4uDh9/4cffpiQkKC/HjVq1MCBA7dt29a+ffu0tDR959dff+3z+QLn0VvHGtEQAgAso4kIlmEKXd4c5EtNTQ3ezMjIyMjICN6TnJxcywrQECLGhDWOCiBMEaxQb/vvIQ/LAABiGj1CAIB1WJgXABDLLJ9r9DxgaBQAENPoEQIALBPBwzJ2j4zSEAIArBVtMUKGRgEAMY0eIQDASrVJkK9XaAgBAJaJJEZod7tJQwgAsEwk6RN29yCJEQIAYho9QgCAdaJwZpna9giffPLJyy+/PCMjo1u3bnPmzAns37hx41VXXdW0adP+/ftv27atbioJAIgOEazKGzUN4Y8//vj8889v3rz5oYceGjVq1OrVq4UQPp9v6NChAwYM2L17d+/evQPrHwIAEC1q2xA+88wzvXv3TktLKyoq6t69+6pVq4QQy5YtO3r06Pjx4xs2bPjnP/95+/bta9eurcvaAgDqNe3nx2XCYXeXMOyHZY4cObJhw4Zf/OIXQojvvvuuc+fODodDCBEXF9exY8ctW7ZYX0cAAOpMeA/LeL3eESNGDB48uHfv3kKII0eONGrUKHA0OTn58OHDNb7R7/d7vd5zqSgAwF4ej8ftdqvLRGMeYRg9Qp/PN3LkSI/H849//EPfk5qaevLkyUCBY8eONW3atMb3aprmcvGEKgBEsZCtYJSqbePk9/vHjh27Z8+eDz/8MD4+Xt+Zk5OzadMmv9+vaZrX6/3+++/btm1bZ1UFANR74SfU2x0irHWPcOzYsZ9//vlLL7104MCBHTt2HDlyRAjRt29fl8v1j3/8w+fzvfDCC5mZmT179qzL2gIA6jct/C+71bYhXL58eWVlZVFRUUFBQUFBweuvvy6EcDqdc+fOnTJlSsOGDd9888233nrL9oWGAQAIS22HRjdu3Fjj/h49eqxfv966+gAAolg0PizDAywAAMtEMum23S0hDSEAwFJ29/DCxeoTAICYRo8QAGCZ+vEcaHhoCAEAlonGGCFDowCAmEaPEABgnShcmJeGEABgHS3soU7bh0ZpCAEAlokkob5ualJ7xAgBADGNHiEAwFLECAEAMUsTmmZ7yxYmhkYBADGNHiEAwDKsPgEAiG11MMfaunXrnn322aNHj95www133HGHucBPP/30+OOPb9q0KS8v74EHHmjcuLEQory8/NNPP123bl3Tpk3vv/9+xfkZGgUAWEaLiOKEZWVl/fr1+8UvfnHPPfc89thjU6dONZe56aabSktLx48fv23btttuu03f+dFHH82aNWvz5s3z5s1T15keIQCg/po6dWp+fv4f/vAHIcRTTz31yCOPjB49OrjAt99+u3LlykOHDiUmJl5++eXNmjX74YcfcnJyRo4cOXLkyJkzZ7744ovqS9AjBABYRo8RhvWlHkpdu3Zt79699de9e/fetGnTmTNnpAJdu3ZNTEwUQjRu3PjSSy/96quvwqozDSEAwDpa+F9KBw4cSE1N1V83bdpUCFFeXn62AnoZqUBIDI0CAGz2ySeftGnTJnjPbbfd9pe//EUI0aBBg0AXUH/RqFGj4JINGjSoqKgIbJ4+fbphw4ZhXZ2GEABgoUjWI/zVr3718ssvB+9MS0vTX1x88cW7du3SX+/cuTMpKUnvFwa0bNkyUEAIsWvXrpYtW4ZVAYZGAQCWiWBkVBOiQYMGrY30FAghxPDhw+fOnXvs2DEhxNSpU2+66Sa9oX333Xc3bNgghLjmmmsOHjz4+eefCyE+/vjjysrK/Pz8sOpMjxAAYB2r1yMsLCzMz8/v1KlTZmbmiRMnPvnkE31/SUlJcXHxZZddlpSUNHny5KKioo4dO27ZsuWVV16Ji4sTQnz22WdDhw6tqqqqrKxMTU295pprZs+eXeMlaAgBAPWXw+GYPn367t27jx07lpub63Q69f2ffvqp2+3WX996662DBg364Ycf2rVrF4gg9unTZ/v27YHzBAqb0RACACwTMkE+Mq1atZL2NGjQIHgzOTm5W7duwXvcbneTJk1qc3IaQgCAZaJxrlEelgEAxDR6hAAAy2iaiCB9oo4qU0s0hAAAS9k91BkuhkYBADGNHiEAwErR1iGkIQQAWCeC9AlihACAC4jVM8ucB8QIAQAxjR4hAMAy0ZhQT0MIALBMHU2xVqcYGgUAxDR6hAAA60ThwzI0hAAAy2jhx/xsH0mlIQQAWEYTmhZmFy/c8pYjRggAiGn0CAEA1okgRmg3GkIAgHXCzyO0veFkaBQAENPoEQIALMOk2wCA2EaMEAAQyzQRQY+wjupSW8QIAQAxjR4hAMAyUTjDGg0hAMBCUdgSMjQKAIhp9AgBAJbRtAjSIUifAABcMFihHgAQ04gRAgAQXegRAgAsE0lCvd1dQhpCAIBlNFafAAAgutAQAgBiGkOjAADLRLQMUx3VpbZoCAEAlokgRmh3O8jQKAAgttEjBABYJwoT6mkIAQDWCT9GaHuQkKFRAEBMo0cIALCMFn4Hz+4OIQ0hAMA6kaRP2B0kPE8NocPhmPnMY0v/NeP8XC7a7d69++KLLw5/Ta/YdebMmePHjzdr1szuikSTQ4cOJSUlNWjQwO6KRJNdu3ZlZ2fbXQvbfP3rXz/66KPqMv2y08I9beKOzF3x8ZFWygKa3+8/P1favXt3dXX1+blWtKusrIy39WMRdfx+v8fjiYuLs7si0aSqqsrtdvP/Vlhi/HfzoosuSkxMtPy0fr+/qqrKxht7/hpCAADqIZ4aBQDENBpCAEBMoyEEAMQ0GkIAQExzTpo0ye46xDS/37927doZM2a8//77paWlnTp1crl+zmn56aefXnzxxXfeecfv9+fk5Nhbz/rp4MGDs2bNSkpKSk9P1/ds3779ueee++ijj9LT0zMzM+2tXn3j9XpnzZo1ffr0NWvWpKen6zfN7/fPnj172rRp27Zt69y5c+DjByGEx+N58803Z86cuWbNmpYtWzZp0kTff+TIEf13UwjB7+YFgB6hzfbv3z9s2LAjR460bNny1Vdf7devn9frFUJ4PJ4rrrji66+/bt269V133fXaa6/ZXdP6aOzYsePHj1+6dKm+uWfPnh49elRVVaWmpl511VWrV6+2t3r1itfrvfbaa19++eUWLVp4PJ6VK1fq+x944IEnn3yybdu28+fPv/HGG+2tZH3zm9/8ZsqUKZdeeumJEyfy8vK2bdsmhKiqqurTp8/69esvueSS0aNHT5s2ze5q4pz5YSuPx+PxePTXJ06cSEhI+Oqrr/x+/9y5czt27FhdXe33++fPn9+6dWv9NQLmzJlzyy23XHnllX//+9/1PRMmTLj11lv11w8//PCwYcPsq129M2XKlLy8vMCHTXfs2LGGDRtu3LjR7/efPn06OTl53bp1NlWw3qmurna73evXr9c3+/bt+8ILL/j9/tmzZ1966aU+n8/v97/33ntt27bVXyN60SO0mcvlCgxGeb1er9fbsGFDIcSyZcv69evncDiEEAUFBTt37ty7d6+dFa1nfvzxx0mTJj333HPBO5ctW1ZYWKi/LigoCPQUIYRYuHDhbbfdtnDhwmeffXbVqlX6znXr1jVu3Dg3N1cIkZiY2KdPH25agMPhaN++/ddffy2EOHbs2Pbt2/UbtWzZsv79++sTERQUFGzbtm3//v021xXnhoawHhk3btzgwYPbtWsnhCgrKwvEveLj45OTk8vKymytXf1y7733jh8/XppTLfimZWRkHD58uKqqyo7a1Uc7d+78xz/+MW/evBMnThQVFU2ZMkUIUV5eHrhjQohmzZrxNz3Yu++++/DDD7dv375Vq1Zjx47t27evMH7MkpKSGjZsyO9mtCMwXl9MmjRp3bp1S5Ys0TddLlfwjHTMHxZswYIFZWVlI0eOlPa7XC49wiqE8Hq9TqfT6XSe78rVVw6Ho2PHjnpAq1u3biNHjrznnnvMH7NYnj9MUl1dPWLEiKKiojFjxuzcufO//uu/evbsmZ+fH/wxE0J4vV5+N6MdDWG9UFJS8q9//euzzz5LTU3V92RlZZWWluqvjx49eurUqebNm9tXwfpl5syZu3bt6t69uxBi69atu3bt2r9//yOPPJKVlRXo0JSWlmZmZtIQBrRo0aJTp07669zc3MOHD+sfqrKyMp/Ppw/Cl5aW9urVy9Zq1iPffvvtV199tWzZMpfL1aFDh1tuuWXatGn5+fnBH7PDhw9XVFTwuxntGBq133PPPTd16tRFixYFD/QNHjz4ww8/PHnypBDi7bff7t69O8kAAU888cTcuXNfeeWVV155pW3btrfeeusdd9whhBg8ePDcuXP9fr8QYu7cuYMHD7a7pvVIUVFRIDT45ZdfZmdnN2jQoHv37m63+9NPPxVC7N+/f+XKlQMHDrS1mvVI06ZNPR7Pvn379M3t27enpaUJIQYPHrxgwYJTp04JId5+++1f/vKX+n5ELybdttkPP/zQtm3b7Ozspk2b6nuefvrpq6++Wghx4403bt++vUuXLgsWLJgzZ06/fv1srWk9ddVVVxUXF48ZM0YIcezYsT59+jRr1iwlJWXlypXLly+P5RVzJKdPn+7bt29iYmLr1q3nz58/depU/R+F6dOnT5gwYeDAgcuWLSsqKnrmmWfsrmk98vvf//69994bNGjQjh07Nm/evGLFiosvvtjv9w8ZMmTPnj15eXkLFiyYO3eu/guL6EVDaLOKiopNmzYF72nTpk1KSooQwufzLV269MCBA3369GnRooVNFazvtm7dmpqaGnh4oaKiYtGiRZWVlf37909OTra3bvWNx+P57LPPTp482atXr+DRvK1bt65bt65169Y9e/a0sXr108aNGzdv3tykSZM+ffoEViAK/G5eccUVWVlZ9tYQ546GEAAQ04gRAgBiGg0hACCm0RACAGIaDSEAIKbREAIAYhoNIQAgptEQAgBiGg0hACCm0RACAGIaDSEAIKbREAIAYtr/A/FTwezU7UcoAAAAAElFTkSuQmCC",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 1
}
],
"cell_type": "code",
"source": [
"using DFTK\n",
"using StaticArrays\n",
"using Plots\n",
"\n",
"# Unit cell. Having one of the lattice vectors as zero means a 2D system\n",
"a = 14\n",
"lattice = a .* [[1 0 0.]; [0 1 0]; [0 0 0]];\n",
"\n",
"# Confining scalar potential\n",
"pot(x, y, z) = ((x - a/2)^2 + (y - a/2)^2);\n",
"\n",
"# Parameters\n",
"Ecut = 50\n",
"n_electrons = 1\n",
"β = 5;\n",
"\n",
"# Collect all the terms, build and run the model\n",
"terms = [Kinetic(; scaling_factor=2),\n",
" ExternalFromReal(X -> pot(X...)),\n",
" Anyonic(1, β)\n",
"]\n",
"model = Model(lattice; n_electrons, terms, spin_polarization=:spinless) # \"spinless electrons\"\n",
"basis = PlaneWaveBasis(model; Ecut, kgrid=(1, 1, 1))\n",
"scfres = direct_minimization(basis, tol=1e-14) # Reduce tol for production\n",
"E = scfres.energies.total\n",
"s = 2\n",
"E11 = π/2 * (2(s+1)/s)^((s+2)/s) * (s/(s+2))^(2(s+1)/s) * E^((s+2)/s) / β\n",
"println(\"e(1,1) / (2π)= \", E11 / (2π))\n",
"heatmap(scfres.ρ[:, :, 1, 1], c=:blues)"
],
"metadata": {},
"execution_count": 1
}
],
"nbformat_minor": 3,
"metadata": {
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.8.4"
},
"kernelspec": {
"name": "julia-1.8",
"display_name": "Julia 1.8.4",
"language": "julia"
}
},
"nbformat": 4
}