{
"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.465892e+01 1.716706e+01\n",
" * time: 0.0030698776245117188\n",
" 1 6.564707e+01 1.018560e+01\n",
" * time: 0.009173870086669922\n",
" 2 5.855772e+01 1.117976e+01\n",
" * time: 0.0223538875579834\n",
" 3 4.332336e+01 1.005994e+01\n",
" * time: 0.04030489921569824\n",
" 4 3.284038e+01 9.111023e+00\n",
" * time: 0.05827593803405762\n",
" 5 2.865867e+01 8.092787e+00\n",
" * time: 0.07396793365478516\n",
" 6 1.425741e+01 2.245998e+00\n",
" * time: 0.08939886093139648\n",
" 7 1.319442e+01 4.157865e+00\n",
" * time: 0.15542101860046387\n",
" 8 1.159360e+01 3.283415e+00\n",
" * time: 0.16910600662231445\n",
" 9 1.061723e+01 3.256131e+00\n",
" * time: 0.18208599090576172\n",
" 10 9.941905e+00 4.446929e+00\n",
" * time: 0.19498491287231445\n",
" 11 8.725013e+00 2.252228e+00\n",
" * time: 0.2078559398651123\n",
" 12 7.578585e+00 1.827392e+00\n",
" * time: 0.22092700004577637\n",
" 13 7.007404e+00 2.221164e+00\n",
" * time: 0.23110604286193848\n",
" 14 6.568958e+00 1.793942e+00\n",
" * time: 0.24155402183532715\n",
" 15 6.366024e+00 1.768311e+00\n",
" * time: 0.2521989345550537\n",
" 16 6.194181e+00 1.054662e+00\n",
" * time: 0.2627260684967041\n",
" 17 6.146899e+00 1.848130e+00\n",
" * time: 0.270676851272583\n",
" 18 5.984643e+00 1.201720e+00\n",
" * time: 0.31085991859436035\n",
" 19 5.847840e+00 1.291968e+00\n",
" * time: 0.32152891159057617\n",
" 20 5.774691e+00 1.191705e+00\n",
" * time: 0.32973599433898926\n",
" 21 5.769751e+00 1.187114e+00\n",
" * time: 0.337831974029541\n",
" 22 5.737326e+00 9.649648e-01\n",
" * time: 0.3481559753417969\n",
" 23 5.732744e+00 1.169599e+00\n",
" * time: 0.3560829162597656\n",
" 24 5.687310e+00 1.049935e+00\n",
" * time: 0.36429500579833984\n",
" 25 5.668944e+00 7.582185e-01\n",
" * time: 0.37228989601135254\n",
" 26 5.645227e+00 6.893850e-01\n",
" * time: 0.38269901275634766\n",
" 27 5.628278e+00 5.089864e-01\n",
" * time: 0.3931128978729248\n",
" 28 5.612760e+00 5.401876e-01\n",
" * time: 0.40117788314819336\n",
" 29 5.600424e+00 4.030253e-01\n",
" * time: 0.4115869998931885\n",
" 30 5.594951e+00 5.280318e-01\n",
" * time: 0.4194190502166748\n",
" 31 5.594137e+00 4.196875e-01\n",
" * time: 0.45024800300598145\n",
" 32 5.584352e+00 2.181642e-01\n",
" * time: 0.45862388610839844\n",
" 33 5.580412e+00 3.504218e-01\n",
" * time: 0.46901798248291016\n",
" 34 5.574590e+00 2.336442e-01\n",
" * time: 0.4771120548248291\n",
" 35 5.571008e+00 2.342273e-01\n",
" * time: 0.48513102531433105\n",
" 36 5.568578e+00 1.322742e-01\n",
" * time: 0.49295902252197266\n",
" 37 5.568511e+00 2.314762e-01\n",
" * time: 0.5009419918060303\n",
" 38 5.567572e+00 1.536611e-01\n",
" * time: 0.5113060474395752\n",
" 39 5.567490e+00 2.450983e-01\n",
" * time: 0.5191559791564941\n",
" 40 5.565619e+00 1.869437e-01\n",
" * time: 0.5271279811859131\n",
" 41 5.564192e+00 1.357510e-01\n",
" * time: 0.5350699424743652\n",
" 42 5.562722e+00 9.982441e-02\n",
" * time: 0.5430610179901123\n",
" 43 5.562696e+00 1.114903e-01\n",
" * time: 0.5508790016174316\n",
" 44 5.562371e+00 1.259793e-01\n",
" * time: 0.5586450099945068\n",
" 45 5.561942e+00 7.301107e-02\n",
" * time: 0.5819690227508545\n",
" 46 5.561517e+00 8.751699e-02\n",
" * time: 0.5902938842773438\n",
" 47 5.561340e+00 9.421214e-02\n",
" * time: 0.5985429286956787\n",
" 48 5.560961e+00 5.570225e-02\n",
" * time: 0.606442928314209\n",
" 49 5.560941e+00 7.654729e-02\n",
" * time: 0.6144099235534668\n",
" 50 5.560832e+00 5.123415e-02\n",
" * time: 0.6247670650482178\n",
" 51 5.560830e+00 7.739870e-02\n",
" * time: 0.6326990127563477\n",
" 52 5.560722e+00 5.386054e-02\n",
" * time: 0.6404340267181396\n",
" 53 5.560654e+00 3.889105e-02\n",
" * time: 0.6485168933868408\n",
" 54 5.560592e+00 4.027462e-02\n",
" * time: 0.6562769412994385\n",
" 55 5.560574e+00 3.126074e-02\n",
" * time: 0.6641368865966797\n",
" 56 5.560548e+00 2.557753e-02\n",
" * time: 0.6744730472564697\n",
" 57 5.560522e+00 1.780308e-02\n",
" * time: 0.6848578453063965\n",
" 58 5.560503e+00 2.154328e-02\n",
" * time: 0.692755937576294\n",
" 59 5.560489e+00 1.565576e-02\n",
" * time: 0.7159719467163086\n",
" 60 5.560476e+00 9.880145e-03\n",
" * time: 0.7241270542144775\n",
" 61 5.560472e+00 1.277672e-02\n",
" * time: 0.7343220710754395\n",
" 62 5.560471e+00 1.186710e-02\n",
" * time: 0.7420790195465088\n",
" 63 5.560469e+00 6.997567e-03\n",
" * time: 0.7525720596313477\n",
" 64 5.560468e+00 7.583083e-03\n",
" * time: 0.7605910301208496\n",
" 65 5.560466e+00 6.826277e-03\n",
" * time: 0.7683570384979248\n",
" 66 5.560465e+00 5.535828e-03\n",
" * time: 0.7787759304046631\n",
" 67 5.560464e+00 4.052430e-03\n",
" * time: 0.7864749431610107\n",
" 68 5.560464e+00 3.954774e-03\n",
" * time: 0.7967419624328613\n",
" 69 5.560464e+00 2.129903e-03\n",
" * time: 0.8073618412017822\n",
" 70 5.560463e+00 1.554788e-03\n",
" * time: 0.8152399063110352\n",
" 71 5.560463e+00 2.222399e-03\n",
" * time: 0.8229219913482666\n",
" 72 5.560463e+00 1.449457e-03\n",
" * time: 0.8435440063476562\n",
" 73 5.560463e+00 1.239992e-03\n",
" * time: 0.8518309593200684\n",
" 74 5.560463e+00 1.009786e-03\n",
" * time: 0.8597838878631592\n",
" 75 5.560463e+00 1.467127e-03\n",
" * time: 0.8677308559417725\n",
" 76 5.560463e+00 1.165432e-03\n",
" * time: 0.8757069110870361\n",
" 77 5.560463e+00 7.346353e-04\n",
" * time: 0.8836190700531006\n",
" 78 5.560463e+00 7.531918e-04\n",
" * time: 0.8915729522705078\n",
" 79 5.560463e+00 1.128620e-03\n",
" * time: 0.8994300365447998\n",
" 80 5.560463e+00 8.181891e-04\n",
" * time: 0.9071450233459473\n",
" 81 5.560463e+00 5.714887e-04\n",
" * time: 0.9175879955291748\n",
" 82 5.560463e+00 6.330130e-04\n",
" * time: 0.9255180358886719\n",
" 83 5.560463e+00 5.163145e-04\n",
" * time: 0.9334168434143066\n",
" 84 5.560463e+00 3.415747e-04\n",
" * time: 0.9413769245147705\n",
" 85 5.560463e+00 3.076041e-04\n",
" * time: 0.9517979621887207\n",
" 86 5.560463e+00 6.823492e-04\n",
" * time: 0.9722809791564941\n",
" 87 5.560463e+00 4.112421e-04\n",
" * time: 0.980604887008667\n",
" 88 5.560463e+00 3.820202e-04\n",
" * time: 0.9887959957122803\n",
" 89 5.560463e+00 3.173862e-04\n",
" * time: 0.9967339038848877\n",
" 90 5.560463e+00 2.687744e-04\n",
" * time: 1.0045499801635742\n",
" 91 5.560463e+00 1.452375e-04\n",
" * time: 1.014922857284546\n",
" 92 5.560463e+00 2.460757e-04\n",
" * time: 1.0253448486328125\n",
" 93 5.560463e+00 1.932301e-04\n",
" * time: 1.035735845565796\n",
" 94 5.560463e+00 1.084828e-04\n",
" * time: 1.046151876449585\n",
" 95 5.560463e+00 9.775272e-05\n",
" * time: 1.059061050415039\n",
"┌ Warning: Negative ρ detected\n",
"│ min_ρ = -1.9499668537205443e-18\n",
"└ @ DFTK ~/work/DFTK.jl/DFTK.jl/src/densities.jl:5\n",
" 96 5.560463e+00 9.775270e-05\n",
" * time: 1.240494966506958\n",
" 97 5.560463e+00 9.775270e-05\n",
" * time: 1.2812418937683105\n",
"e(1,1) / (2π)= 1.7391793934279238\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3deXwURd748eo5ciIJCQmRgEQIcrkcCwi7oCKQKJdGQDSru4DoI7jHw+4+PKyuootnVp/1QlZXERQVWVTwxaKuJyAol3gAEYzchIRD5CbJTGZ+f7Q7v+nqUJMeOnSG+bxf+WO6u6a70pmkUvXtb5UWDAYFAADxyuV0BQAAcBINIQAgrtEQAgDiGg0hACCu0RACAOIaDSEAIK7REAIA4hoNIQAgrtEQAgDiGg0hACCuec7OZQKBwF33Tp/65z+fncvFOr/f7/GcpR/NOaO2ttbtdjtdi1hSW1vrcrk0TXO6IrEkzn830xK9Ects2XLo1Cm/pdPu2PHdvn0bbrvttmjrdabO0k+0urq60qd9uOPg2bkcAMB2IzucH7HM9dcv/uqrAxZPvLGo6IiDDSFDowCAuBa/fXwAgO00TcTccDs9QgBAXKNHCACwjaZpVp/ACgYd7kLSEAIAbMXQKAAAMYQeIQDAPjH4sAwNIQDANlHECIXTMUKGRgEAcY2GEAAQ1xgaBQDYJoqEeqdHRmkIAQD2iSJG6PjDNQyNAgDiGj1CAIB9tNhLqKchBADYRrM+1On40CgNIQDANtHkETrdEhIjBADENXqEAAD7RBEjjFT+0KFDL7/88uHDh0eMGNGjRw9zgUAgsGDBgtLS0q5du44cOVLvkr711lv79u0LlcnOzi4qKqrz/PQIAQB20qx/KRw7dqx3796rV68WQgwcOPD99983l5k4ceJf//rX1NTU6dOn//73v9d3btmy5fP/ePDBB998883TXYIeIQCg8XrllVfOP//8l19+WdO0zMzMBx98sKCgILzA7t27X3rppR07duTk5Nxwww0dOnS48847s7Oz//d//1cvUFVV1bJly5tvvvl0l6BHCACwjab9+LyMBco+4UcffTRkyBB9tHPIkCHLly/3+XzhBZYuXdq1a9ecnBwhxAUXXJCfn79ixYrwAm+++WZ6evrll19+ukvQEAIA7BPFwKhybHTv3r3Z2dn66xYtWgQCgcrKyvAClZWVLVq0CG3m5OTs3bs3vMCsWbMmTJigeJaVoVEAgG30DqG1t2ja2rVrx4wZE77zqquu0gcz3W53IBDQd+ov3G53eEmXyxUqIISora0NL7Bjx47ly5fPmTNHUQEaQgCAw1q2bHndddeF7+nYsWPoUEVFhf567969Ho8n1EE0FxBCVFRUtGzZMrT5wgsvFBYWtm7dWnF1GkIAgG2iWH1CCJGbmys1hCFDhgx5/PHH7777brfbvWjRosLCQo/HI4T45ptvMjIyWrRoMWjQoAkTJmzbtq1t27bffPPN7t27Q+HAQCDw4osv/u1vf1NfnRghAMA+dscIx4wZ43K5rrzyyttuu62kpGTatGn6/nHjxs2bN08IkZ2d/fvf/76goGDy5MlDhgy544470tPT9TLvvffeiRMnhg8frq4yPUIAQOOVlJS0YsWKJUuWHD169K677goNcj799NOhZ2QeeOCBoUOHbtq0qbi4uE+fPqH35ufnf/jhh4mJiepL0BACAGyjiQjpEHW9JYKkpKRRo0ZJO3v16hW+2a9fv379+kll8vPz61MBGkIAgH2sxwidnnObhhAAYKMYXI+Qh2UAAHGNHiEAwDZRJNQ7PjZKQwgAsE0UeYROt4MMjQIA4hs9QgCAfWLwYRkaQgCAbX5cWMniWxqoMvVEQwgAsE0MdgiJEQIA4hs9QgCAfaLoEjrdhaQhBADYJ4oYYQPVpN4YGgUAxDV6hAAA22jWE+SdfmiUhhAAYB9Ns54O4XRLSEMIALBPDD4sQ4wQABDX6BECAGwTgx1CGkIAgH2imWLN6aaQoVEAQFyjRwgAsE8Mjo3SEAIAbBNFHqHjDSFDowCAuEaPEABgnxica5SGEABgH2KEAIB4pmmxN9coMUIAQFyjRwgAsI0mNKsJ8jGTUH/ixIm77rqrf//+AwYMeOSRR2pra/X95eXlv/rVr3r16jVhwoT9+/c3WD0BALFAs/7ltPo2hH/605+WL1/+9NNPl5SUvPDCC0888YS+/9prr23atOmcOXOEEMXFxQ1USwAAGkh9h0bXrVt36623duvWTQgxZsyYdevWCSFWr15dVla2cuVKr9f75JNPZmVllZaWdu7cuQHrCwBoxKJ4WMbxTmF9e4RFRUXz5s3bsWPHN99889ZbbxUVFQkhvv766x49eni9XiFEampqly5dvvrqqwasLACgcdN+nHbbGmfrXN+G8LbbbquqqurZs2efPn3y8/OvueYaIcT+/fvT09NDZZo1a3a6MGEwGAyFFQEAscjn80UudA7HCEePHt23b98DBw4cPHjQ4/H85je/EUI0bdr05MmToTLHjx9PS0ur8+2aprnd7jOvLgDAKfr437mnXg1hMBhctmzZL37xC5fLlZCQcMMNNyxdulQIkZeXV1ZWppcJBALbtm3Ly8trsKoCABq9aMZGHa5yvRpCTdM6dOjw9ttvCyECgcDbb7/dqVMnIURhYeGxY8eWLFkihFiwYEFiYuKll17aoNUFADRmMTgyWu+nRp9//vlf/vKXc+bMqa6uzs7Ofu2114QQiYmJs2fPHj9+fHp6+vHjx19++WXGPwEAsaW+DWHfvn3LysoOHDjg9XrDH5AZNmxYeXl5ZWXl+eef7/EwTw0AxLcGmHR73759s2fPPnbs2IgRI/r27Wsu4Pf7586du2XLli5dutx4440ulyu0f/78+Rs2bGjWrNmoUaPy8/PrPL+1uUazsrLCW0Gd1+tt3bo1rSAAIKrkCVVLePTo0UsuuaSsrCwtLW3YsGHvvPOOucyECROeffbZ3NzcJ554Qn+WUwhRU1NTUFAwc+bM9PT0H3744dNPPz3dJWi9AAC2iWKFenX5l156KS8vb9asWUKItLS0Bx98cMiQIeEFduzY8c9//nPXrl1ZWVkjR47Mz8+fNm1aTk7OzJkzjx8//tlnn0Xsp7H6BACg8Vq2bNmVV16pvy4sLFy5cmVNTU14gU8++aRr165ZWVlCiNzc3IsuumjlypVCiCVLlowdO/add9556qmnvvjiC8UlaAgBAPaJ4rFRZY+woqIiOztbf92iRYtgMFhZWXm6AnqZvXv3CiG2b98+Y8aMN998c9++fVdeeaU+J3adGBoFANgmmsRATVu7du2YMWPC9xUUFNx6661CCI/HE5qYzO/3C1Nev8fjCQQCoU2fz6cX0DStW7dus2fPFkJ069btj3/847hx4+q8Pg0hAMBhLVu2vO6668L3dOnSJXRI7+EJIcrLyz0eT3j/Ty9QXl4e2ty7d2/Lli2FELm5uaGTXHzxxeXl5X6/v854IQ0hAMA2Uaw+oQmRm5srNYQhw4cPLykpufvuuz0ezxtvvDF06FA9Yf3LL7/MysrKzc0dPHjwzTffXFZW1r59+w0bNlRUVAwYMEAIcc0113z00Uf6SdasWZOfn3+6p2ZoCAEAjdeoUaNmzJhxxRVXtGvX7l//+te///1vff9tt91WXFw8efLk5s2b/+lPfxo8ePCQIUOWLFkybdq0pk2bCiEmTJgwZ86coUOHtm7deuHChS+++OLpLkFDCACwj/UYobp4YmLi0qVLP/jggx9++OHhhx/OycnR9z/33HP6k6JCiGnTpg0dOrS0tPS2227r0aOHvrNp06Zr1qx57733Tp06dc899+jjpXWiIQQANGoJCQlDhw6Vdnbt2jV8s1evXr169ZLKJCYmjhgxIuL5aQgBALaJxRXqaQgBALZyumGzioR6AEBco0cIALBNFAn1jq/MS0MIALBNNHmETg+l0hACAGyjiQjLKtX9JkcRIwQAxDV6hAAA+0SxQr3TaAgBAPYhRgg0dkHFVqMm/61w+m8HcM6gIQQA2Eazng5BjxAAcA4hRgg4wNpop/L4Gb3Zmsh/K4wl5EsHVSdgHBVOicU8QtInAABxjR4hAMA2UUyx5niXkB4hACCu0SNEDAhGCM2pgoSW3hr5eKSq1F/Q9F+wvB20cFg6KAcQlQFF4fx/5ICTaAgBALaJ5mGZhqlJ/dEQAgBso4nYW4aJGCEAIK7RI0QjEDmqZyUKKB1VHpaPmq4dtBJyjBDUk2jyueTFa6QtdUhReTFNqpqpcNAYRCQNEdGLIqHe6Q8YDSEAwD7WY4SOoyEEANgmijxCYoQAADiJHiEcECEZTx3WM20Hg6qwn7WjAWXFTOWtzT0qx/zM/wUbUwOl8i5jGE+TCktHpUtHildKoVOpvBRhjLWBL0CNhhAAYJuohkYbqC71xdAoACCu0SNEw4gwgmgxh0E9vBkwbAeUo52BgHRUtSmEkHcoa6IWcbjSpR4LNW66lCOl8lHjf7x1/AMu71IN0pJrAYUoZpZx/DNDQwgAsJXTDZtVNIQAANtomuV0CNInAABwEj1C2CTCtGdBxVF1zM+8J6DerLVyVFnYvEeOKZ5BjFCK+QlTYE/edCs31YWVZxbmIKKpgOGoMtdCnjku1kbJcKaiWH3C6Q8JDSEAwD5RzDXqNIZGAQBxjR4hAMA2mtDkdVTq8yZH0RAiWpaCglamPZPidsIUqKutNeQG1vqNR/1WjhpPZb50ba0yRhhQxT7lvEFlXqAQwi1F8jyGTbfbMH7j9kibmoWjbnkoyBVUbcojR1JQ0FLIUDj+Rw8NK5oV6p3+SNAQAgBsFWMdQhpCAEDjtmfPnueff/7w4cNFRUUDBgwwF6ipqZk1a1ZpaWnXrl3Hjx/v8XiEEKtXr166dGmozMSJE9PS0uo8Pw/LAABso0VBecLDhw/36dPn4MGD+fn511133VtvvWUuM3bs2Ndee61r166zZ8+eOHGivnPp0qXz5s374T8CgdMuLkOPEPWjni1UWAsKWsrtE6bAnt8Y9vP7jEd9tYbNGsNRn1xYCijKvyrqIKI0kWlQmqJTSEFBQ2Ep1U9EigJ6vIZNr3HTk6Aq7PG6DfX0yHdYulZQSkMUUhqi9G4LIUNBouG5zva5Rl988cWOHTvOmDFDCJGcnFxSUnLNNdeEF9i6deuiRYvKy8szMjJGjBiRl5c3ffr0li1bCiEuueSShx9+OOL16RECABqvTz75ZPDgwfrrQYMGrVq1qqamJrzAypUru3fvnpGRIYTIycnp2LHjp59+qh/avHnz/fffP3v27MOHDysuQUMIALCPFtXX6VVUVGRlZemvs7Ozg8FgRUXF6QroZfQCLVq06Ny5c21t7WuvvdapU6edO3ee7hIMjaJe6phKzL6xUCnDwW8an5QGMH3VhsFPX41hs0Y6Wi0dNQ6NGt8rDZyKSMkYpvQJ9dColOFgGhr1qAc/DcObCYnGwonSUcOmN8F4/41HhRDSWKknaOFfZEsjpcKcX1H/KyEmRLEwr9BWrlxZUFAQvvPqq6/+7W9/K4Twer1+v1/f6fP5hBCJiYnhJRMSEkIF9DIJCQlCiHHjxo0bN07fOXr06JKSkpkzZ9ZZARpCAIBtophhTROiffv2U6dODd/Ztm1b/UVubm55ebn+ury83Ov1hvf/hBAtW7bcs2dPaHPPnj25ubnSJfr06fPxxx+frgIMjQIAHJadnT3YKNQQXnPNNQsXLtTjgvPnzx8+fLjb7RZCrFmzZteuXUKIwsLCrVu3lpaWCiHWr1+/f/9+PcXixIkT+hl8Pt+SJUsuvvji012dHiEAwDb1SIio4z2Kg9dee+3f//73/v375+XlLV269P3339f3//a3vy0uLp48eXKzZs2mTZtWUFAwePDg99577/7772/SpIkQ4vLLL2/WrFlmZua6deuaNWt2xx13nO4SNISom2mBoQjbZxQU9EkxPzlQJ4X9aqoMm9VVfsPmKWlTeq/hqBxQNF1anV8hL9ukTAxwyTFCeTxGTpAwZkRIYb+EJMMvb2Ky4ahfPhphtalgwBg1lH/8ckxRQQ4ZmtMnpPyKoDFkSMzwHGBr+oTX6/3ggw+WLVt29OjRZ555Rn86VAjx0ksvhV5PnTp1+PDhpaWld9xxR8eOHfWd//73vz///PNjx4794Q9/6N27t6J5piEEADRqbrd74MCB0s4OHTqEb3bp0qVLly7hezIzMwsLC+tzfhpCAIBtopl0u2FqUn80hAAA22iasJw+4XRLSEOI/6gjVTDsoHmKNSk6ZgyuWQoKSrl9UhhPmMJ+VdLmCcPmqZPGoyd9xlOp4otS0qEwRQ0bNEYoBQWl1MBEOexnqHlStTd805+iXm1K/sWXEyKD6r8MqpChvBZdpL9x8h9BJmCLddHlTziK9AkAQFyjRwgAsE0MdghpCAEANrI7j/AsoCHEj9TLKtWxDFPAQuKgNH2oFHiTgoJSCFAIceqEIc536rihwMnjNcajxsInpJChKumwpo4YoWGPFOwMGGNv0j2Tfrtd8ipL8i+/1ytlChpjhMmG39akFMOmFN30+wwhw9pIC10pA8RC0zzStnHLZdw0nMxl+hOnybtUyzY53ldAPLAQI/zyyy+Lioo6dux4xRVXrFy5Ut+5Y8eOUaNGderU6YYbbti7d2/DVBIAEDvsW3ri7KhvQ7h58+aCgoIBAwYsWrRo2rRpKSkp+v6ioqJ27dq9/fbbzZs3v/766xusngCAGKDnEVr9clZ9h0anT59+4403Tp48WQgRmsDm008/3b1794MPPujxeB555JHmzZtv2LDhJz/5SUNVFrZST6KmnkFNRBoaldeUN05UJo1ASuOT0kCoEOLkMcOeE0cNY6EnjEdPHjMcPWkcKZVyLaRhWGnyNmHK9DANjUpZB8ZlmIy/39KS9OahUY9yaDTJODRaXWUcGq0xpk9YqaeZVDNpMFPalEY65cKmv3FB+ezKBe6ZgC3WRDHXqOWYot3q2yNcu3Zthw4dxo8fP2TIkCeeeKK2tlYIsXHjxu7du3s8HiFEcnJyly5dNmzY0ICVBQDAbvXtEe7evfuvf/3rjBkz0tLSJk6c+MMPP9x7770HDhxIS0sLlWnWrNn+/fvrfHswGNTbTgBAjPL5fF6vN0KhxhH2s6S+PcK0tLRbbrll2LBh/fv3/8tf/vLKK68IIdLT00MLPgkhjh071qxZszrfrmmavoIUACBGRW4Fo4oROt5w1rdH2K5du1DnLz09/fjx40KICy+8cMuWLfrOQCDw3XffhZZSROyRo4LGLVNQSQ4KGhMJpBiVT7mOkhSok7IjhCkoePyIavOEFCM0RhDVM66ZY4Q+ZYzQNDOZKkaoyVOsWUufMC0gZQwKystFqepppg5nSokf0lxx0qbmNtRESq4Q5iBihNQN5VHADvXtEY4dO/bVV189efJkbW3tCy+8oK+IMXjw4Orq6jfeeEMIMXfu3CZNmvTr168BKwsAaNxi8KHRevcIJ0yYsHbt2jZt2ni93m7dus2ZM0cIkZCQ8PLLL//qV7/63e9+5/V6X331VZe8KCcAIJ7E4Bxr9W0IPR7P888//9RTT2malpSUFNo/aNCg3bt3HzlyJD093fFHYAEAsMraFGvJycnmnS6X63TPyKBxOYPEQXOQyZQ4qFpoSZqoTF5W6YRqyjRhyhQ0xQirDYWPqvIIpRnXTIE3cx6htICRKnsywhRrLlUcTgjh8RoGVBKqjTHCasNvq6lixsnelEFB0yRnphxHeVM1OZy0KccXXXJNXNJNs5JWKJfmf+/GyHoeodM/SOYaBQDYJpqgn9P/0NAQAgDsE4MxQp5tAQDENXqEcSTC0jvKxMGAOUZYK8UIrUwuaozMnTop5RGa5xpVZQrKU48eNU49elyVRyglDkrrQwnTAlKmUJyhsHSPTDN2Gjal5DwhhMdj2OOvUd1SefpQueLSpSNNc2rcI8166klwqTa9UgTRUBW3KRQaMEYNLaUVEiJs/KKZa9TpnyQNIQDANjE4MsrQKAAgvtEjjF+mfAlpU5UYICKmTxiH9XzyukvGKdZOqjIchGl4U5o1TT4q50tIc6qp8iX8pqFRaeo403pGhsJB402VBnxMC9abrmUcLK2tNWyqUzUk6mFYaTBTCOFNMIyFehMNNy0h0Tj3W6Lxphnf6/EaM0y8piH5CJ801ein410HRBaDXUIaQgCAbbQo1hekIQQAnDOiyCN0uh0kRggAiG/0CM9pERImVEfPMH3CJ6dPGDeN6ROmkKEcI5TmYFPHFE0LLfkVm1K+hPRdCCFq5SwFdaBOlUAh/Zus1ZryCtzK+dvUs6bJAUgpHcI4eVuivDhoYpLhtiQmG/4y1CRLgVXD0YQkw03zygkn8n/b0vflkudYUyehMONao0eMEAAQ16LII3R6wQaGRgEAcY0eIQDANky6jcZFHRJUz6kmB8MixghrVfOB+Y3LMKlnXJPCeEKIKjnOZzybtJRSlepaPuW8ZeYYoTp7L6jO5pPSCqUFhky3NBgwFIiUKWj4vqQ1njxe6Y4ZgoJVSfIdTpLusLw6lXEFKONPU76HysxLEemjFVTPuGb8i6kOJwL1REMIALCNFoPrERIjBADENXqEAADbRJNQ7/SgNg1hPLEQz4oQMhSm4FnAONeoOq3Qp5yJVArymfdYiwIqMwVr5aQ307epnnPVSh5hUJNChkJmuriivL/GsO1zS4mbyjsW6Q77lLfUpwysSp8EcwZqhFxMK59Sp0fUUIdolmFyuiWkIQQANGrbt29/5plnjhw5UlRUdNVVV5kLVFVVPf3001u2bOnSpcukSZMSEhLCj86fP//gwYO//vWvT3d+YoQAAPto1r+UDh069LOf/SwQCPTu3Xvs2LFvvPGGucxNN9307rvvXn755QsXLrzlllvCD3322We/+c1vHnzwQcUl6BHGL9OIlGpMyjy9l7wMk7wqk2r0zK8eKfXJA3d+4x6pgJyq4Vc9zS9NmSYP8JrnTJPzKSyshSSnT0jH5KnF5KXbjSOjcs1dtarvWr7D8h2T77BUwFKSiXyHAxGGRpVTxcmfw6A0TZ3qrWgc7M4jnDNnTrdu3R555BEhhNvtLikpGTVqVHiBb7/9dsmSJZWVlWlpaYWFha1bt37ggQdat24thKiurv71r389bdq0hx9+WHEJeoQAANtoPz4uY+lLZeXKlQMHDtRfDxw4cN26ddXV1eEFPvvssx49eqSlpQkhsrKyOnXqtGrVKv3QvffeO3r06E6dOqnrTEMIAGi8Kioqmjdvrr/OysoKBoOVlZXhBSorK0MFhBDZ2dkVFRVCiC+//PK9996bMmVKxEswNAoAsE9Uq0+sXLmyoKAgfN+wYcMmT54shEhMTKypqdF36i+SkpLCSyYmJvp8/3/Nmerq6qSkJL/ff8stt8yYMcPr9Ua8Pg1hHDGFv4KKw/LkYXUEz1QTZckRRGVyhXqzrrdLs7sZL62c4kv97H5dMT9LQUEV6b2aPHuYHDVUZ7AErHzX6htYVwH1D0hVk4iT8wk5I0W5spJUWCNk2Nhp1vMCNSHat28/derU8J3t27fXX7Rq1WrPnj366927dyckJGRlZYWXzM3N3b17d2hzz549ubm533777ddff33TTTcJIU6dOrV///527dotXbpUjx1KaAgBAA7Lzs4ePHhwnYeKioruuuuuu+++Oykpad68eVdffbXL5RJCLF++vFWrVm3bti0oKBg/fvzXX3/dtWvX1atXf//991dccYXH49m8ebN+hpUrV/7P//zP+++/n5OTU+claAgBALbRNOsJ8sry11xzzXPPPXfJJZfk5eWtXbv2ww8/1PdPmTKluLh48uTJ6enpDzzwQGFh4WWXXbZs2bKSkpKUlBQhRNu2bfWS27Zt83g8oU0zGkIAgH3sXqHe4/G8/fbba9asOXTo0Msvv9y0aVN9//z58/UnRYUQ//3f/z1s2LDNmzeXlJRceOGF0hn69ev36aefqi5hsb6IG5GCZ5GWbVIdlbP3lNEv854I6/hIl1ZPHRcxFHoGQUE185k1ZTwyUgTRsCmvHhX5DhtDjMofkKWfdcRPTqTVwhBjoplrNFIBl8vVt29faWdeXl74Zn5+fn5+fp1vT05ObtOmjer8kSoAAMC5jB4hAMA2sTjpNj1CAEBco0eIH51hpEY9c6k8kakySlRXVEmVeabMhzTPoRq7ISnVLZW/SWUoVI6MRozbyXfYQk2sit0fD2IXDSEAwE72Trp9FtAQAgBsE02MsIGqUm80hPjRGX4WNXlTU2yaSxu2TFWRfq9MmxZqYl4NSXnUXMBBqlsqf5PyLVXdQHP5CD8gKzWxyvG/iThTducRngU8LAMAiGv0CAEAtonF9AkaQgCAbaKYWcZxNIQ4jchxO+OmyxgocqmOutSbbnnE3uWWTq56u3xp5Tci/SsaNC2NJEcU7YsY1vHHwlIoVL7/hk3TPZHusHxt6Z6rf0CWftYRPzmOx4cAYoQAgLhGjxAAYJtYjBHSIwQAxDV6hHFEDnfJUSZD+Ev+H82ceaaOQhk33R5p06XclK9lersq4iWHu9yGVYK0WmO1A9J3LV1ZXvxI01Szu6mZb6G6gDoKG+m7Vt//iHfYwg9IHcGVNoWIEApVZ0QST4wJtq7LezbQEAIA7BODCfU0hAAA20QRI3QcMUIAQFyjRxi/Is3JaQh/mQM9chRQDgqqgkwer2HTm2Dc9Lqla3mMe6QC0tk8xmvVGoNnAbcxyBcwfmOmmF/AJQUFpZBV9EFC8z/NLuMu+Q5LYT+36ruW77B0xxLkO+w17pF+ItLZ5B+uW9pUfReirs9SOHtnLsXZF0VCveMdSBpCAIB9YjBGyNAoACCu0SOMJ1YWIFI/uy/Ms3Apn7/3qsdCEw3jcglJ8sCdtEferDZs+msM+RK1tcah0YA6/yEglDtMK7/LxRXUU6aJOnJOlMPLCapbmqC8pdJR8x7pJ2Iau1ZVTM7rMH1y1DO0RegcON11QERaHXGWiG9h0m0AwLlCiyKPsGFqUn80hAAA+0QXI3R0AWxihACAuEaP8FymjglaCgrWESNUP80vpTQkqIs8Y/QAAB2mSURBVGJUiUmGz2FisvyxTDLuqanyGzarDUf9PkNYTw4KBlT/eZqHdGprVW+XsinUcdeIYVc5QcJzBlFA+ZYajiaZ7rB0z6WfiHQt6acpJ1coJ3sTdazTpAydRsi1QOMTRUK9Zi0RyXY0hAAA20STR9gwNak/hkYBAHGNHiEAwDZRrUfYQHWpLxrCc5pp4SVpW1FaWm+ojomy1AEtrxTQMm7KESxDzC8pRf5YVlcZY4TV3vBNOShYa2GlJM1VazhVjenbrLVw8qDxlsqzhRnPbQ6eyXFWZVBQiuolpRjuSXKq6mhSqikKm6KKKSbIIUNV8FKdVihMn6WIq1Opjjr9BxTnBhpCAIBtoogROv4PDTFCAEBco0cIALCNJqzHCJ3uEtIQxi91aMaU3SWH2kwTYxo25Zkw5fiWMa3NGPPzVcsTfvpqjJE8OSho2AyapgsNJwfqjN+Fzy2/2e+3cC3l7K3ypJoutzwe41EnDhozBaWwX0oT5eZ5qqNCiORUYxDReHLp5yX9NE2LNEWca9S4qcwjdHrMDNad86tP+P3+X/ziFxMnTgzt+fbbb6+66qoLL7xwxIgRO3futLt6AAA0LGsNYUlJyYYNG1avXq1vBoPBa6+9tk+fPqtWrerYseOYMWMaoIYAgJihPyxj9ctZFoZGN2/evGDBgj/84Q9PPvmkvmfFihX79++fNm2a2+2+7777srKyvvzyy+7duzdMVXGmIiRTSIWVQ4giUvqER06fMK6UZHz0X1o4ye8zDIQKIfz+hPBN9bRnEtOApHEI1yvN1ma6tDI3Q72ok3pONbcpr0B906QcBilBQhrtTG1quGOp5xk2U4ybQojkJsaTy9kUyhnXvIbNyOkTluZUY6Q05kQxxZrTP9j6NoSBQODWW299+umnw8c/N23a1L17d7fbLYRISkrq3LmzvqdBagoAaPwaIEb47bffPv3004cPH7722muLiorMBU6cOPHEE0+UlpZ27dr1d7/7XVJSkhBizZo1b7zxRnl5eUZGxujRoy+77LLTnb++Q6N/+9vfunfv3q9fv/Cd33///XnnnRfaTE9PP3DgQJ1vDwaDfr+/zkMAgJjgyJ/xAwcO9OvXr2nTpoWFhbfffvu8efPMZYqLi1euXFlUVPTBBx+MHz9e37lhw4bMzMzhw4e3bNly+PDh//rXv053iXr1CHft2vX4448vWbJk27Zt+/fvr66u3rZt24UXXtisWbPjx4+Hih09ejQjI6POM2ia5vHwhCoAxLD6/Bm3vUM4Z86c3r1733fffUKIQCDw6KOPFhcXhxfYvHnzBx98sG/fvvPOO2/gwIEtW7bcuXNnmzZtJkyYECqzZ8+ef/3rX8OHD6/zEvVqnPbu3ZuYmDhy5EghxPHjxw8fPlxQULBx48a2bdtu3rw5GAxqmlZbW1tWVtauXbv6nBDOUAYJ5XVQpECOaewgQvqEMW7kTTCcPDHJsFlrTFGQQoDCPLFZhKCg8fF95fpQ1acM9aypkmOEPmPAstavXpXJOMWafA9Vd0wI4fWqVlaSc05SVAkSUlCwSZoUMpTTJ1KaGApIc7BJqzJ51cswRU6fUC7DJJUmSBhrNNPPNPJblMU/++yzAQMG6K8HDBgwduzYqqoqffBTt2rVqh49eujDkxkZGV26dFm9enWbNm1CBXbv3r1ixYo//vGPp7tEvYZG+/btu/U/Hnvssc6dO2/dujU5OXnw4MGBQOCVV14RQjz33HMZGRk///nP63NCAADqo7KyMjMzU3/dvHnzYDBYWVl5ugJ6mYqKCv31vHnzNE274IILfvKTn9x0002nu4TlKdaaNGmSm5urv/Z4PK+99tq0adPS09Mff/zxV1991frDQgCAc4gWzdfSpUt7GZWUlOjnS0pKqqmp0V9XV1cLIVJSUsIvmJyc7PP5QptVVVWhAsXFxcFgsKys7LvvvrvnnntOV2XLcburr7766quvDm32799/27ZtUkcVABCfossL7N69+6OPPhq+J9TjatWq1a5du/TXu3btSkpKat68uVQyVEAIsXv37latWoUXyM/PnzRp0hNPPDF9+vQ6r27PAyy0gucCK2mFwhwjNIbigh5jvp0x80xKvwvUehRHhSn2JtdNmZ+nTs6rSjLmEZpihH5ljNC0KpMqRmjKvDTnEapihEnSukvKPEIpU1AKCkpZhsKURyhdS6qJPMWaR/VJqCNGaCVxEDEnqvUItfT09J49e9Z5dNSoUVOmTLn77rtTUlJefvnloqIil8slhHj//ffbtGlz0UUXFRYW3nzzzZ9//nnPnj1XrFhx5MgRPaa4f//+7OxsIUQgEPjwww8vuuii01WAJzkBAI3X8OHDX3jhhZ/+9KetWrXavHnzhx9+qO+/6667iouLL7rooqZNm5aUlFx11VV9+vRZvXr1//3f/yUnJwshhg4dGgwGs7Ozt2zZct555y1evPh0l6AhBADYx+78CbfbvWjRoi+//PLo0aO9e/cOxf8WLVqUmpqqv540adLw4cO3bNnSuXPnli1b6jtXrVq1cePGgwcPtmzZslOnTop+Kg0hAMBGmtVllSKW1zStR48e0s7zzz8/fLN169atW7cO3+PxeOo50xkNYfyS/j0KWkkrFOY5PI3F5alHjUeDAbdx0xhpq7u+4XVRxd5MixkZs/GMQcGkUxFihOo8QnlVJuU9U8/OKiLnERrDeCmqqUeT5VWZEoxH5V98aRmmSJOLSomDqqCgVkeM0ELiIM+h4yygIQQA2EaLwX9faAgBAPaJwYV5aQjxI9P8a/IOubw0k5lU2jjQ6gm6pMPGrQifQ6kqpswN5exuiYZLS4N+1VWGoVGfaRkmn3GJKGk2ODnTQ7lEvWlSOtPQqHFJeilLQZrnTD3jmrSOkpRrIWVHCPNYqDynmnGgW7nQkrRpHho15UuoR0oRY6JLn2igytST5ZllAAA4l9AjBADYJoqZZRwfBqAhBADYhxghYpgyJhj5XzxlqoCJW3FM0+SPpXplJTlVI0E1p1pNsjFfotpwLSkiKITw+85ejFBKS/Aqv5EEZcjQFFBUTZlmPrkcFDSGXeWFlqSgoLRukjlEqP5oOP03EXGIhhAAYJsoHpZx/L8fGkIAgG1icGSUp0YBAPGNHiHqpp6ATUSag81l7V8sY8iqjuncpJOr1jMyrbtkSA2UgoL+GsNRn0+OEUpBQWmKNdPkcMZlmKTZwpTVFqaoodcrxT6NMcJEVdKhKeanWkdJmG6alCmoDgrKc6pFzgtkErVzXaz9TGkIAQD2icGEehpCAIBtoskjdLoHSYwQABDX6BGiXsz/sQXlUJAqh85SyFDT5NLSdJWWYoR+KbSWZIj5qdMEzXukxMGgsXikGKGhsMs0CadbjsxJyXzqCKKqsJQI6IkUnlRPHxohKKicSlTEXvwIFkXx2KjTaAgBALbRNOsxP6fHRhkaBQDENXqEqB/zRFnS5GJnMFJqGluLsOST5jaMSJpWfjcc9XgNZ/NK6RDGJeYDfvnStfJYqDp9wlhtaTPS0u1uKUvEOIDpVk4sJ40PRzjqNg0+KxNUtAg/LyFOv13Hv/qxNm4GS6JYmNfpDiENIQDAPprQNIv/7DjdDtIQAgBsFINzrBEjBADENXqEiJa02o6lkKE8gZthy/zfmWllH+Oz/i7D+93GkGHAGCMM1LqMm8ajATlGaMqXMG7KxVXk78IUI1QH6iylNEQ4arq0On4p11x+szIo6PQ/+zjbrCfUO46GEABgmyiWYdLMD8idXQyNAgDiGj1CAIB9YvBhGRpC2ORMQobSmUwRrKAyuuYyxu0CLlUYL0LMzxwjlHYEbYsRmgMpLuW3qV6LSjqZKREwwspH8liWHJSVi5++rPN/1OCsaCbdNuXgnmU0hAAA20SRR+j4f0/ECAEAcY0eIQDAVjHWIaQhRANRZp4Fg+qQlSleIIUY5QCksawxHhY0xvFcxtJW8wKls1mLbMiBtwi//eq8Q1NupYWjdcwcq762hYOIdyzMCwBAjKFHCACwTTQJ9Q1UlXqjIYQD5F8TZa6FiJRuEWHgNKi6WFBKWTBdJ2iqnKIi6mWYZHWMT55BDoOcm6G8VKS1kciIQPTIIwQAxLUo5hp1uiEkRggAiGv0CAEAttHq8US06S0OdwlpCNEIRHi23xTnk7fUEUTVdsRJ74PS+y1NqqamDvrVUd7KQavndnpsClDYuHHjk08+eeTIkaKiouLiYnOBo0eP/vWvf92yZUuXLl2mTJmSmpoqhFi9evXChQu3b9+enZ19yy23dOvW7XTnZ2gUANB47du377LLLsvLyxszZszUqVNfeuklc5nrr7++tLR0/Pjx69atGzt2rL7zvvvuS0pKGjNmTEZGxs9+9rOvvvrqdJegRwgAsE106xEqjr7wwgv9+/e/8847hRA1NTUlJSW/+tWvwgts2rRp+fLl+/fvT01N/fnPf56Tk7Nt27a2bdsuXrxYP/OoUaO+/PLLhQsXnq5TSI8QAGAbfWYZS1/qkfs1a9Zceuml+utLL73066+/rqqqkgr06NFDHw5NT0+/+OKL161bJ8La12AwuGfPnvPPP/90l6BHiBgQ4f9Lada0CBHHoOqweXI3+VK2BdMin+gMYoYkAsIxUeQRKlVWVmZmZuqvmzdvHgwGKysr8/LyQgX27dsXKiCEyMzMrKioCD/DjBkzjh07JvUjw9EQAgActnTp0l69eoXvGTlypD4cmpKSEuoCnjp1Sgihd/5CkpOTq6urQ5tVVVXhBebPn//QQw99/PHHycnJp7s6DSEAwDZRxQhF9+7dH3300fCdF1xwgf6idevWO3fu1F/v3LkzOTk5vP8nFdDLtGrVSn+9cOHCyZMn//vf/+7QoYOiAsQIAQD2sR4j1DSRnp7e0ygrK0s/3+jRo19//fXjx48LIebMmTNq1CiXyyWEWLJkSWlpqRCisLCwoqJi1apVQoiPP/74xIkTV1xxhRDi3XffnThx4uLFi7t27aquMj1CxL5IaYgGUkBReapIlzqrCPshPg0dOnTu3Lldu3bNzc3ds2fPhx9+qO+fPn16cXFx586dmzRp8thjjw0fPrxHjx5ffPHFU089lZiYKIS4/fbbq6urr7/+er38mDFjHnrooTovQUMIAGi8XC7X/Pnzt2zZcujQoZ49eyYkJOj733nnnaSkJP31+PHjhw0bVlZW1qFDh+bNm+s7v/jii0AgEDqP3jrWiYYQAGAbTUSxDFPk8uYgX0ZGRvhmdnZ2dnZ2+J60tLR6VoCGEHHG0jgqAIuiWKHe8d9DHpYBAMQ1eoQAAPuwMC8AIJ7ZPtfoWcDQKAAgrtEjBADYJoqHZZweGaUhBADYK9ZihAyNAgDiGj1CAICd6pMg36jQEAIAbBNNjNDpdpOGEABgm2jSJ5zuQRIjBADENXqEAAD7xODMMvXtET700EO9evXKzs7u2bPn/PnzQ/s3btx4+eWXZ2ZmDh48uKysrGEqCQCIDVGsyhszDeH333//xBNPlJaW3n333ePHj1+zZo0QIhAIjBo1aujQoTt37uzXr19o/UMAAGJFfRvCRx99tF+/fs2bNy8qKurdu/fq1auFEMuXLz98+PCUKVOaNGny5z//eevWrevWrWvI2gIAGjXtx8dlrHC6S2j5YZlDhw5t2LDhpz/9qRDim2++6dq1q8vlEkIkJCR06tRp8+bN9tcRAIAGY+1hGb/fP3bs2BEjRvTr108IcejQofPOOy90NC0t7eDBg3W+MRgM+v3+M6koAMBZPp/P6/Wqy8RiHqGFHmEgEBg3bpzP5/vHP/6h78nIyDh+/HiowJEjRzIzM+t8r6ZpHg9PqAJADIvYCsao+jZOwWBw0qRJu3bteueddxITE/Wd+fn5mzZtCgaDmqb5/f5vv/22ffv2DVZVAECjZz2h3ukQYb17hJMmTfrkk0+efvrpffv2bdu27dChQ0KIgQMHejyef/zjH4FA4Mknn8zJyenTp09D1hYA0Lhp1r+cVt+GcMWKFdXV1UVFRQUFBQUFBS+88IIQwu12L1iwYObMmU2aNHn11Vdfe+01xxcaBgDAkvoOjW7cuLHO/ZdccslXX31lX30AADEsFh+W4QEWAIBtopl02+mWkIYQAGArp3t4VrH6BAAgrtEjBADYpnE8B2oNDSEAwDaxGCNkaBQAENfoEQIA7BODC/PSEAIA7KNZHup0fGiUhhAAYJtoEuobpib1R4wQABDX6BECAGxFjBAAELc0oWmOt2wWMTQKAIhr9AgBALZh9QkAQHxrgDnW1q9f/9hjjx0+fPjaa6+9+eabzQV++OGHBx54YNOmTd26dbvzzjubNm0qhKisrPzoo4/Wr1+fmZl5xx13KM7P0CgAwDZaVBQnrKioGDRo0E9/+tPbb7/9/vvvnzVrlrnMddddV15ePmXKlLKysptuuknf+e67777yyiulpaULFy5U15keIQCg8Zo1a9aAAQN+//vfCyEefvjh6dOnT5gwIbzA119/vWrVqgMHDiQnJ/fq1atFixbfffddfn7+uHHjxo0bN3fu3Keeekp9CXqEAADb6DFCS1/qodR169b169dPf92vX79NmzadOnVKKtCjR4/k5GQhRNOmTS+++OLPP//cUp1pCAEA9tGsfynt27cvIyNDf52ZmSmEqKysPF0BvYxUICKGRgEADnv//ffbtWsXvuemm276y1/+IoRITU0NdQH1F+edd154ydTU1KqqqtDmyZMnmzRpYunqNIQAABtFsx7hz3/+82eeeSZ8Z/PmzfUXF1xwwY4dO/TX27dvT0lJ0fuFIa1btw4VEELs2LGjdevWlirA0CgAwDZRjIxqQqSmprY10lMghBBjxoxZsGDBkSNHhBCzZs267rrr9Ib2zTff3LBhgxDiyiuv3L9//yeffCKEeO+996qrqwcMGGCpzvQIAQD2sXs9wsLCwgEDBnTu3DknJ+fYsWPvv/++vr+kpKS4uPgnP/lJSkrKjBkzioqKOnXqtHnz5meffTYhIUEI8fHHH48aNaqmpqa6ujojI+PKK6+cN29enZegIQQANF4ul2vOnDk7d+48cuRIly5d3G63vv+jjz7yer366xtvvHH48OHffffdRRddFIog9u/ff+vWraHzhAqb0RACAGwTMUE+Om3atJH2pKamhm+mpaX17NkzfI/X623WrFl9Tk5DCACwTSzONcrDMgCAuEaPEABgG00TUaRPNFBl6omGEABgK6eHOq1iaBQAENfoEQIA7BRrHUIaQgCAfaJInyBGCAA4h9g9s8xZQIwQABDX6BECAGwTiwn1NIQAANs00BRrDYqhUQBAXKNHCACwTww+LENDCACwjWY95uf4SCoNIQDANprQNItdPKvlbUeMEAAQ1+gRAgDsE0WM0Gk0hAAA+1jPI3S84WRoFAAQ1+gRAgBsw6TbAID4RowQABDPNBFFj7CB6lJfxAgBAHGNHiEAwDYxOMMaDSEAwEYx2BIyNAoAiGv0CAEAttG0KNIhSJ8AAJwzWKEeABDXiBECABBb6BECAGwTTUK9011CGkIAgG00Vp8AACC20BACAOIaQ6MAANtEtQxTA9WlvmgIAQC2iSJG6HQ7yNAoACC+0SMEANgnBhPqaQgBAPaxHiN0PEjI0CgAIK7RIwQA2Eaz3sFzukNIQwgAsE806RNOBwnPUkPocrnmPnr/sn++dHYuF+t27tx5wQUXWF/TK36dOnXq6NGjLVq0cLoiseTAgQMpKSmpqalOVySW7NixIy8vz+laOOaLX/zivvvuU5cZlNfc6mmTt+XsSEyMtlI20ILB4Nm50s6dO2tra8/OtWJddXV1oqMfi5gTDAZ9Pl9CQoLTFYklNTU1Xq+X/7csifPfzfPPPz85Odn20waDwZqaGgdv7NlrCAEAaIR4ahQAENdoCAEAcY2GEAAQ12gIAQBxzX3vvfc6XYe4FgwG161b99JLL7311lvl5eWdO3f2eH7Mafnhhx+eeuqpN954IxgM5ufnO1vPxmn//v2vvPJKSkpKVlaWvmfr1q2PP/74u+++m5WVlZOT42z1Ghu/3//KK6/MmTNn7dq1WVlZ+k0LBoPz5s2bPXt2WVlZ165dQx8/CCF8Pt+rr746d+7ctWvXtm7dulmzZvr+Q4cO6b+bQgh+N88B9Agdtnfv3tGjRx86dKh169bPPffcoEGD/H6/EMLn81166aVffPFF27Ztb7311ueff97pmjZGkyZNmjJlyrJly/TNXbt2XXLJJTU1NRkZGZdffvmaNWucrV6j4vf7r7rqqmeeeaZVq1Y+n2/VqlX6/jvvvPOhhx5q37794sWLR44c6WwlG5tf/vKXM2fOvPjii48dO9atW7eysjIhRE1NTf/+/b/66qsLL7xwwoQJs2fPdrqaOGNBOMrn8/l8Pv31sWPHkpKSPv/882AwuGDBgk6dOtXW1gaDwcWLF7dt21Z/jZD58+ffcMMNl1122d///nd9z9SpU2+88Ub99T333DN69GjnatfozJw5s1u3bqEPm+7IkSNNmjTZuHFjMBg8efJkWlra+vXrHapgo1NbW+v1er/66it9c+DAgU8++WQwGJw3b97FF18cCASCweCiRYvat2+vv0bsokfoMI/HExqM8vv9fr+/SZMmQojly5cPGjTI5XIJIQoKCrZv3757924nK9rIfP/99/fee+/jjz8evnP58uWFhYX664KCglBPEUKIt99++6abbnr77bcfe+yx1atX6zvXr1/ftGnTLl26CCGSk5P79+/PTQtxuVwdOnT44osvhBBHjhzZunWrfqOWL18+ePBgfSKCgoKCsrKyvXv3OlxXnBkawkZk8uTJI0aMuOiii4QQFRUVobhXYmJiWlpaRUWFo7VrXH73u99NmTJFmlMt/KZlZ2cfPHiwpqbGido1Rtu3b//HP/6xcOHCY8eOFRUVzZw5UwhRWVkZumNCiBYtWvA3Pdybb755zz33dOjQoU2bNpMmTRo4cKAwfsxSUlKaNGnC72asIzDeWNx7773r169funSpvunxeMJnpGP+sHBLliypqKgYN26ctN/j8egRViGE3+93u91ut/tsV66xcrlcnTp10gNaPXv2HDdu3O23327+mMXz/GGS2trasWPHFhUVTZw4cfv27f/1X//Vp0+fAQMGhH/MhBB+v5/fzVhHQ9golJSU/POf//z4448zMjL0Pbm5ueXl5frrw4cPnzhxomXLls5VsHGZO3fujh07evfuLYTYsmXLjh079u7dO3369Nzc3FCHpry8PCcnh4YwpFWrVp07d9Zfd+nS5eDBg/qHqqKiIhAI6IPw5eXlffv2dbSajcjXX3/9+eefL1++3OPxdOzY8YYbbpg9e/aAAQPCP2YHDx6sqqridzPWMTTqvMcff3zWrFkffPBB+EDfiBEj3nnnnePHjwshXn/99d69e5MMEPLggw8uWLDg2WefffbZZ9u3b3/jjTfefPPNQogRI0YsWLAgGAwKIRYsWDBixAina9qIFBUVhUKDn332WV5eXmpqau/evb1e70cffSSE2Lt376pVq4YNG+ZoNRuRzMxMn8+3Z88efXPr1q3NmzcXQowYMWLJkiUnTpwQQrz++us/+9nP9P2IXUy67bDvvvuuffv2eXl5mZmZ+p5HHnnkiiuuEEKMHDly69at3bt3X7Jkyfz58wcNGuRoTRupyy+/vLi4eOLEiUKII0eO9O/fv0WLFunp6atWrVqxYkU8r5gjOXny5MCBA5OTk9u2bbt48eJZs2bp/yjMmTNn6tSpw4YNW758eVFR0aOPPup0TRuR3/72t4sWLRo+fPi2bdtKS0tXrlx5wQUXBIPBa665ZteuXd26dVuyZMmCBQv0X1jELhpCh1VVVW3atCl8T7t27dLT04UQgUBg2bJl+/bt69+/f6tWrRyqYGO3ZcuWjIyM0MMLVVVVH3zwQXV19eDBg9PS0pytW2Pj8/k+/vjj48eP9+3bN3w0b8uWLevXr2/btm2fPn0crF7jtHHjxtLS0mbNmvXv3z+0AlHod/PSSy/Nzc11toY4czSEAIC4RowQABDXaAgBAHGNhhAAENdoCAEAcY2GEAAQ12gIAQBxjYYQABDXaAgBAHGNhhAAENdoCAEAcY2GEAAQ1/4fPULB14CTYrMAAAAASUVORK5CYII=",
"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.3"
},
"kernelspec": {
"name": "julia-1.8",
"display_name": "Julia 1.8.3",
"language": "julia"
}
},
"nbformat": 4
}