{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Numerical calculations for \"Origins of Cosmological Temperature\"\n",
"\n",
"This SageMath notebook performs numerical calculations for the paper *Origins of Cosmological Temperature* and the supplemental note *Calculations for \"Origins of Cosmological Temperature\"*.\n",
"\n",
"It makes graphs of the two anharmonic potentials and does some arithmetic.\n",
"\n",
"Section headings are as in the supplemental note.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3.5 Action in dimensionless variables"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAaEAAAHWCAYAAADejza7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xd4VNXaPv57UoAkJPQakK6oiC0qilgQpQhYKR4rx4IN8diP2H0Vjx55fRU5iKjYEBRQEJGmUkQQRWnKQVFAkZLQUkggbX5/PN/9W3tmT5KZycysXe7PdXFlZRZ7ZrHZs5+9us/v9/tBRESkQZLuAhARkXcxCBERkTYMQkREpA2DEBERacMgRERE2jAIERGRNgxCRESkDYMQERFpwyBERETaMAgREZE2DEJERKQNgxAREWnDIERERNowCBERkTYMQkREpA2DEBERacMgRERE2jAIERGRNgxCRESkDYMQERFpwyBERETaMAgREZE2DEJERKQNgxAREWnDIERERNowCBERkTYMQkREpA2DEBERacMgRERE2jAIERGRNgxCRESkDYMQERFpwyBERETaMAgREZE2DEJERKQNgxAREWnDIERERNowCBERkTYMQkREpE2K7gLU5I03gP/7P+D334FDh4C0NKBrVyDp/4VPvx8oLgZ++w0oKwPq1wcKC/WWmYiIwuPz+/1+3YUIx/33A//+NzBhAnDbbdb84mLglluANWuATZsSXz4iIoqcY5rjli+Xn4MGhc5PTwfGjgU6dkxcmYiIqHYcEYQKC6WGc8IJQJs2gXkFBSqdnAx06JDYshERUfQcEYSWLwfKy4H+/QNf37YNGDNG/Z6WBgwbltCiERFRLTgiCH31lfw0B6GDB4HRo4Fzz1WvNWoE9OqV2LIREVH0HDEwIScHWLsWOOUUGRW3b5+MhgOA3FygaVO95SMioujYviaUnw/8+KPUcFavBlatAn79FXj/faB795oD0PjxwJIl1tfHjZP+pTZtrPlbt6q8UaOsxw4erPKDTZmi8mbNCswrLFR5V19tPXbECJW/b19g3ty5Km/SJOuxHTtKXr9+1rz771fHbt4cmLdypcobO9Z6bE6O5OXkWPPGjlXHrlwZmLd5s8q7/37rsf36SV6ogSSTJqlj584NzNu3T+WNGGE99uqrVX7wUP1Zs1TelCnWY428wYOteaNGqfytWwPzlixReePGWY/t3l3yzj7bmvfkk+rYNWsC8zZsUHnmZmdD796S17WrNW/8eHXswoWBebt2qbyRI63HDhmi8o8cCcybOlXlTZ0amHfkiMobMsT6viNHqvxduwLzFi5UeePHW4/t2lXyeve25o0Zo47dsCEwb80alffkk9Zjzz5b8rp3t+Z57R4R6tpNBNvPE1q2DKisBHr2DHz9uOOsF+QvvwBHH61+/89/gLp1gXnzZNCCuamuoAD46y9JB3/RKipU3oED1jLl5an8YIcOqbzi4sA8v1/l7d1rPXbfPpVfWRmYV1Ki8oqKrMf+9RdQWgo0b27NO3BAHVteHph35IjKMw/yMOzeXfW/tbpzWF5e/TnMzZX8OnWseUVF6tiSksC8ykqVF/wlBOS8GvnBdfziYpV36JD1WCOvbVtrnvkcVlQE5tV0DnftknLVq2fNy89Xx5aWBuaVlam8gwetx+7ZI/mZmda8wkJ17OHDgXnm63v/fuux1V3f5nMYfH0DKi8vz5q3f3/V5/DwYZUXap7fzp3yeoMG1ryDB9WxZWWBeaWlKi8/33qscX0HX7+A9+4RGRmhyxtvtg9CRn/QWWcFvn7kSODTVnEx8PTTwLvvyu+HDklAuuAC+X327MDjs7KA7GxJ160bmJecrPIaNbKWqVkzlR8sI0PlpacH5vl8Ki9UDa5JE5WfFFRHTUtTefXrW4/Nzq76AmvUSB2bEvQ/XreuysvKsh7bsmXgT7PqzmFKSvXnsHlzyQ8VhOrXV8empQXmJSWpvCZNrMc2baryfb7AvPR0lRfqC2fkNWtmzTOfw+TkwLyazmGrVvJ3Qp3DBg3UscHnIjVV5TVsaD22RQu5sYa6HjIz1bHBwc98fTdubD22uuvbfA6Dr2+g+nPYuHHV57BePZUXKqi2bi031hYtrHkNG6pjU1MD8+rUUXmhAljLlhIAQ50Hr90jQl27iWD7PqGTTwbWrZMngFD/2YZJk+QCveSSxJWNiIhqx9ZBaP9+eaI49lhg48aq/15ZmdR4vvxSPekvXQrMnAkcf7xUnd97D/jmm8SUm4iIwmPr5rilS0P3BwV7/HHpYDQC0PLlsrTPqlVSxRwzxto/QERE+tk6CM2bJz/PPDN0/m+/AS+8IE1xv/4qr/n9Mgpn1CjVxpmfD5x3XtyLS0REEbJdc9y+fcBFF8nP7dvltTZtAjskDx+WkUHG6JE+fYBFiyS9ejVwxhnAn3+q4ZHdugEvvgj07Zu4fwcREdXMdjWhJk2s8yUisXWrBCwjAO3dK0O3e/aUPqHgUXZERKSP7SerRuq44wKHLr74oryWkcGBCUREdmO75rhYePJJmS9Qrx5w2mkys/+MM4DrruMq20REduLKIERERM7guuY4IiJyDgYhIiLShkGIiIi0YRAiIiJtGISIiEgbBiEiItKGQYiIiLRhECIiIm0YhIiISBsGISIi0oZBiIiItGEQIiIibRiEiIhIGwYhIiLShkGIiIi0YRAiIiJtGISIiEgb1wWhzz4D3n1XdymIiCgcrglCfj/wn/8A48YBFRW6S0NEROFI0V2AWPH5gNtuA/bs0V0SIiIKl2tqQkRE5DwMQkREpA2DEBERacMgRERE2rhmYAIATJ4MLFoENGwIpKcDQ4fqLhEREVXH5/f7/boLQURE3sTmOPr/+f2cY0XuVVkJlJfrLgUF82wQGjcOaNNG/ixZEpi3davKGzXKeuzgwSo/2JQpKm/WrMC8wkKVd/XV1mNHjFD5+/YF5s2dq/ImTbIe27Gj5PXrZ827/3517ObNgXlffgk0aADUrQukpcnP444DbrgB+O47ICdHjsvJsb7v2LHqfVeuDMzbvFnl3X+/9dh+/SSvY0dr3qRJ6ti5cwPz9u1TeSNGWI+9+mqVX1gYmDdrlsqbMsV6rJE3eLA1b9Qolb91a2DekiUqb9w467Hdu0ve2Wdb8558Uh27Zk1g3oYNKm/MGOuxvXtLXteu1rzx49WxCxcG5u3apfJGjrQeO2SIyj9yJDBv6lSVN3VqYN6RIypvyBDr+44cqfJ37QrMW7hQ5Y0fbz22a1fJ693bmjdmjDp2w4bAvG++keb5OnWAlBSgXj35/7jpJuCHH+T/pE0beS2Y1+4Roa7dRHBNn5DPV3VeqAbHggLgr78kHfxFq6hQeQcOWI/Ny1P5wQ4dUnnFxdZyGHl791qP3bdP5VdWBuaVlKi8oiLrsX/9BZSWAs2bW/MOHFDHmp8EZ8+Wi7qgIPDvb9okf95+WwJTSYn1PYHqz2F5efXnMDdX8uvUseYVFaljgz+7slLlBX8JATmvRn7w/3txsco7dMh6rJHXtq01z3wOg2uLR46ovOBzCcgNd+9euQEGy89Xx5aWBuaVlam8gwetx+7ZI/mZmda8wkJ17OHDgXnm63v/fuux1V3f5nMYfH0DKi8vz5q3f3/V5/DwYZUX/PAAADt3yusNGljzDh5Ux5aVyU+/H3jzTeCf/5RzbKiokEC1YQPwxhtAVpb8nwVfv4D37hEZGaHLG2+uCUKR9mxlZQHZ2ZKuWzcwLzlZ5TVqZD22WTOVHywjQ+Wlpwfm+Xwqr2lT67FNmqj8pKA6alqayqtf33psdnbVF1ijRurYlP/3P/7KK8Do0YHnrUED4KijpBZj3BBLSuR8NGxofd/qzmFKSvXnsHlzyQ8VhOrXV8empQXmJSWpvCZNrMc2barygx9M0tNVXqgvnJHXrJk1z3wOk5MD8+rWVXlZWdZjW7WSv9OypTWvQQN1bPC5SE1VeaHOf4sWcoMNdT1kZqpjg4Of+fpu3Nh6bHXXt/kcBl/fQPXnsHHjqs9hvXoqL1RQbd1abqwtWljzGjZUx6amygPQLbcAb70V+PeaN5fjf/5ZBcGCAjnvod7Xa/eIUNduInBgggdNmADccYf6ffBg4F//Us06hw8D77wDPPywqm00aQLMnx+6WY7ILoqLgWHDAptxBw8Gnn5aNbkVF8tI2n/+U9VE+vQBPv00dG2V4stVQWjmTHkSysuTJ5FQbfte9+OPQI8eqqYzZox8QUM1Z+7eDVx8sbSdA/Kk9PXXwAknJK68ROEqL5fv/Oefy+916gDvvRe6fwqQPrgLLlDNdUOHAtOmVd+0T7HnmoEJmzcD27bJRXjjjXIxUaDCQvmiGQFo9Gjgf/6n6i9dy5bAV18B55wjvxcUAIMGhW4DJ9LtwQdVAMrMlJp7VQEIAE49Vf6+0TT74YfSj0SJ5ZogVL8+8OqrMqKrsDB0u6fX3XYbsGWLpHNygOefr/mYrCz5op52mvy+fTtw882R98ERxdP8+Wp0V2qqNK2df37Nx515ZuD+Y/ffH3pAAMWPa4JQdjbwxBNAz57AgAHS3kvKl18C778v6cxMqSmGGhQQSno6MGOG6oCdORN47bX4lJMoUnv2ANdfr35/4QXg3HPDP/6yy4Brr5X0gQPAo4/GtnxUPdcEocOHZVjxsmVSvR48OPQQUi+qqADuuUf9/tJLQKdOkb3HUUcFNlXcfTfwyy+xKR9RtPx+mWaQmyu/9+8P3HVX5O/zr3+pEWWvvSZ9p5QYrglC77wjNaAePaRqfvbZ8pOAjz8G1q2T9CmnyETUaFx6KXDnnZI+ckQm6bFZjnT64APVD9SihUwEjWZgQatWqgbk9wMPPBCzIlINXBOEUlICJ7mddhrQrp2+8tiF3w88+6z6/ZlnrPMLIvHcc2oy58KF1hnfRIly6FBgsJg8uXZ9waNHAx06SHrxYuD772tXPgqPa4LQddcBGzdKk9Fbb0nV+tRTdZdKvwULVNPCqacCffvW7v0yMqQ5z3D33aFXHyCKt+eeUysEXHwxMHBg7d6vbl0ZYWcYO7Z270fhcdU8IbLq0wf44gtJz5wJXH557d/T75emT6O587HHZA00okTZtk0mVx85IqPhNm4Ejj669u97+LCssbZrlzTr/f470L597d+XquaamhBZbd+uAlCnTtKnEws+H/Dyy2oJoP/9Xw5rpcR64AG1ntvo0bEJQICsmGCsJuL3W5f+odhjEHKxt99W6REjatcXFKxLF1mJGJC+uBdeiN17E1Vn3Trgo48k3bx57IdUm78rb77J7U3ijUHIpSor1XYFPp/0mcXamDFqYcdXXpFlfoji7amnVHrMmNgvvNm6tfQxAcCOHdatMCi2GIRcavlyte9Nnz6htyeorTZtZBUGQFbbZkcuxdu6dWpEZqtWsnpHPNx4o0pPnhyfzyDBIORS5k3bQm3+FisPPaSWo5840bpZGVEsmWtBDz1k3eojVgYMUFtvzJkjqzJQfDAIuVBRkWozb9AgdgMSQmnRQnXklpbK+n1E8bBpU2JqQYCMuDMmdZeXq+8TxR6DkAvNmKHm7gwfHr+nRcPo0Wqk3MSJXC6J4sO8/fR998X/uv7b31R6xoz4fpaXMQi5kPkLY17YMV6ys2UjMUA2wTOvSkwUC3v2yNJcgAxEMEZmxlO3bmro97JlbJKLFwYhlzl0SJYcAWSUzxlnJOZz//EPlX7pJRmdRxQrr76q9sG65ZbEbEXt8wFXXilpv1/WYKTYYxBymcWL1SS+QYNiOzeoOqeeqja/++9/uXgsxU5xsWxJD0izbzSrZEfLCEIA10mMFwYhl5kzR6UTvb25ebuI//3fxH42uddHH0kzLyDNvvGYblCVk05SCyEvXSqDfii2GIRcpKJCdpQEZNh0796J/fyBA4HOnSW9eDGwfn1iP5/cybyBojESM1F8PhmuDUhzoNHUTbHDIOQiq1cDeXmS7ttX1sFKpORkGSlnePnlxH4+uc+GDcDKlZLu1k32C0s0Y/UEAPjss8R/vtsxCLmIsbkXIP1BOtxwg+o0njYtcI8nokhNmqTSI0dGt2FdbZ1/vnqgmzePGznGGoOQi3z1lUpfeKGeMtSvr+ZXHDokO18SRaO4WA33T0sDrrlGTznS0yUQAcDOnTJplmKHQcgliouBb7+VdJcusq6bLuaZ7OYnWaJIfPghkJ8v6WHDgIYN9ZXlggtU2tgehWKDQcglVqwAysokbTy16XLKKWpX2zVr1M6uRJEIborTiUEofhiEXMLcFKc7CAGBtaHXX9dXDnIm84CEE05I3KTrqnTvDjRtKumvvpL15Cg2GIRc4ssvVfq887QV4/931VVqde3331dr2RGFww4DEsySktSUh4ICqeFTbDAIuUBhIfD995I+9li1BL1OWVkSiAD50nIBSApX8ICEq6/WWx4Dm+Tig0HIBZYvV1sQJ3qCanXMG4MZi08S1WTmTDUgYfhwvQMSzMzN3CtW6CuH2zAIuYDd+oMMPXrISD1AyvjHH3rLQ85gfmAxP8jo1rmz6hdatYqL9MYKg5ALmIPQuefqK0cwnw+47jpJ+/3Ae+/pLQ/Z344dqqmrUyfgrLP0lsfM5wPOPFPS+/cDv/yitzxuwSDkcPn5agi0eQSPXZgnGL7zDmebU/Xef19dI9deq39AQjBzUDRG71HtMAg53OrVqlmgVy+9ZQmlfXtVO9u8GfjuO63FIRvz+wOb4q69Vl9ZqmLUhAAGoVhhEHK4VatU2vwFsROjSQ4A3n5bXznI3n78Efj5Z0n37Al07Ki3PKHk5MhCvQDwzTd6y+IWDEIOZw5COlYYDseVV8pQW0CWYuFEPwrl/fdV2o61IADIyABOPFHSP/+sRvFR9BiEHMzvV0GoaVN7PjkCMmdo4EBJ790bOLGWCJAm5Q8/lHRKCjBkiN7yVMfoF/L71XqNFD0GIQfbskVG6QBSC7JbJ67Z8OEqPW2avnKQPa1cKSPjAFkBvnFjveWpjrnZm01ytccg5GBOaIoz9O8PZGZKetYs4MgRveUhe5k+XaXNDyx2xMEJscUg5GBOCkJpacCll0o6Px+YP19vecg+KiqAjz6SdJ06wCWX6C1PTdq3V0tjcdJq7TEIOZjRHu3zAaedprcs4TDWkgPYJEfK8uXA7t2S7t8faNBAb3lqYp60WlCgRvRRdBiEHKq0VJa7B4BjjlFbattZnz6qrX/OHK6sTcLcFDdsmL5yRMLcJLd6tb5yuAGDkENt3CiBCFAbyNldaqoM1wZkpeS5c/WWh/QrL1crrKelAYMG6S1PuMzfubVr9ZXDDRiEHMrYugGQCXROwVFyZPbllzJsH5Bh/PXr6y1PuE46SaW5c3DtMAg5lHlTLafUhADgnHOAVq0kPW8eJ/t5nROb4gBpVj7qKEmvXcvBCbXBIORQxtOXzxf4VGZ3ycnA0KGSLi0FPv5Yb3lIn9JSGa4PSA1owAC95YmU8b0rKgJ+/11vWZyMQciBKiqkTwiQ5e6N+TdOYX7iNYbmkvcsWgQcPCjpwYPV0k5OcfLJKs0muegxCDnQli1ASYmkjXWsnOSMM4DsbEkvXswmOa9yalOcwRyEODghegxCDrRunUo7MQglJQGXXy7p0lLgs8/0locS7/Bh4JNPJN2gAdC3r97yRIODE2KDQciBnB6EAOCKK1R65kx95SA95s8HCgslfdllQN26essTjaOOAho1kjSDUPQYhBxo/XqV7t5dXzlq4+yzgebNJf3555y46jXmASnGQBWn8flUk9zu3WrVB4oMg5ADGTWhBg2Adu30liVayclqLbmSEq4l5yXl5cCnn0o6MxPo3VtveWrD3CTHfqHoMAg5zP79wJ9/Srp7d3tv31ATNsl50/LlwIEDkh4wwJlNcQYOTqg9BiGHMdaLA5zbH2Q4/3zVpj53Lrd38IrZs1Xa7itm14RBqPYYhBzGPCjBqf1BhtRUmR8CSCf14sV6y0Px5/erUXGpqc6boBrs6KPl3wFwNe1oMQg5jBtGxpmxSc5b1q8Htm+X9Pnn23/bhpqkpkogAoDNm6W/iyLDIOQwRnOczwd066a3LLFw4YVq0crZs4GyMr3lofgyakGA85viDMcdJz9LS2UiOUWGQchB/H5V5e/YEUhP11ueWKhXD7j4Yknv3w8sXaq3PBRf5iBkNMU63fHHqzSb5CLHIOQgf/6p5tMce6zessQSm+S8Yft21XmfkwO0aaO3PLFiDkI//aSvHE7FIOQgmzaptNEE4Ab9+0uNCJBJjFwW353Mo+KMOWJuYP4uMghFjkHIQcxVfTcFofr11dphe/YA336rtzwUH+amODcFoS5dOEKuNhiEHMR8gbupOQ4IvCmZn5jJHfbvB5Ytk3SnTu56iOIIudphEHIQc3Oc24LQwIGyujYQ+MRM7vDZZ7IPFiAPHE5e6SMUjpCLHoOQQ5hHxrVp47yN7GrStKksagrI0+R//6u3PBRbbm2KM3CEXPQYhBxizx613pabmjLMzPNG2CTnHiUlwIIFkm7WDDjzTL3liQeOkIseg5BDuHVknBmDkDt98YWaWjBokKyg7jYcIRc9BiGHcPOgBEOnTmoViFWruD+LW7hpwdKqdOmiguvmzXrL4jQMQg7h1uHZwYz+Ar9f7TlDzlVRAcyZI+n0dFmmyY1SU2UVEwD45Re5fik8DEIOYe6od2tNCGCTnNusWgXk5kq6b18gLU1veeLJGKZdXAzs3Km3LE7CIOQQv/4qPxs3Bpo00VuWeDr1VCA7W9KLFwNFRXrLQ7XjhaY4wzHHqPQvv+grh9MwCDlASYnaTbVLF71liTefT92sjhxRo6rIefx+WYYJkDlgAwfqLU+8GTUhgP1CkWAQcgDz5De3ByEg8ImZE1eda9Mmde326uXuGjwQGIRYEwofg5ADGE1xgDeC0HnnAVlZkv7sM+4x5FRuXbC0KgxC0WEQcgBzEDJf6G5Vp47a9vnAAWD5cr3loei4cQO76rRuDWRkSJrNceFjEHIAr9WEAC5o6nQ7dwKrV0u6e3egQwe95UkEn089JG7dKuvIUc0YhBzAi0Gof3+1PP4nn3DehdMYc4MAbzTFGYwgVFEhgYhqxiDkAEYQat5c9ZW4XVYW0Lu3pP/4A1i3Tm95KDJuX7C0KuwXihyDkM0VFQG7dknaK7UgA0fJOVNBAfDll5I+6ijgpJP0lieROEw7cgxCNue1QQlmgwerNPuFnOPzz9WIxksucd/eQdUxT1hlEAoPg5DNebE/yJCdDZx2mqTXrgW2b9dbHgqP10bFmZm/o7/9pq8cTsIgZHNeDkJA4E3M3NlN9lRaCsybJ+mGDYFzztFbnkRr2BBo1EjSv/+utyxOwSBkcwxCKs1+IftbskT6hABZpscY4eglnTrJzz//5DDtcDAI2Zw5CHXurK8cuhx/vFoif+lStbss2ZOXm+IMRhCqrAS2bdNaFEdgELI5IwhlZ6vZ2F5iXtC0okI19ZD9VFaqJtO6dWXrBi8yHpoANsmFg0HIxgoLgbw8SRtPV17E1ROcYc0a4K+/JN2nD5CZqbc8upi/qxycUDMGIRszV+W9sOxJVc46S63A/PnnssUD2Q+b4gRrQpFhELIx87IfXg5CKSlqL5qiIjURkuzFqKX6fMCgQXrLohNrQpFhELIxBiGF237b25YtwE8/SbpHD6BlS73l0Sk7W1aCBxiEwsEgZGMMQspFFwH16kl6zhzpBCf78NreQdVJTgbat5f0779z8d2aMAjZmDkIGRe1V2VkABdeKOldu4Dvv9dbHgrk1QVLq2I0yRUXA3v26C2L3TEI2ZgRhFJTZcMsr+PEVXvKzQVWrJB0167eW+MwFPYLhY9ByKb8fjU6rl07qeJ73cCBajFM9gvZx9y5qsmJtSDBEXLhYxCyqf37ZZ4QwP4gQ4sWwJlnSvrnn6UznPRjU5wVa0LhYxCyKQ5KCI2j5Ozl0CFg0SJJt2qlVj33OtaEwscgZFMclBBadasnvPUWMH16YsvjdQsXAocPS3rwYCCJdxQA0oRuCN6CJC8PuPlmTro28JKxgS1brBcka0KhHX20dH4D0hluLGsEAN99Bzz0EIfExsvGjdYbKldJCC0zE2jcWNLB52zrVmDyZOC//018ueyIQUizbdvkxjpsWODrDEJVM252lZXSKW4YMkTO5+rVWorlaitXAieeGLgSQnm5Ov/16wO9e+spm10ZLRg7dsi5Cn6dK2wLBiGN/H7g2mvl5/LlgU/wXDeualX1C51zjszUnzYt8WVysxUrJMBUVgIlJer1r7+WATQAMGCArJxNitEkV1GhFnYFgGbNgLQ0BiEDg5BGr70mX2RAvsyrVqk8oyaUni4XLSlnnCEj5QDpkygulnRyMjB0qPQLVVToK5+bLFsmK2IbzcX79qk8NsVVz9yXaw44Pp/kMQgJBiFNtm0D7rlH/Z6cDLz3nqTNm2G1b6/mxpBISpJOcECezBcvVnnDh8uKCkZwp+h99ZXsCVRaqmrpBw7I1AG/XwWhlBSpCVGg6gYnMAgpDEIa+P3AiBFAWZl6raIC+OADeW33bvXkyaa40KpaPaFHD/nys0mudhYvBvr3lwAUvE7f1q3AunXqxnr++UDDhokvo91VVRMC5BplEBIMQhq89hqwZElgZyUgT5lffMFBCeG44AK10+zcuar5zeeTQR4zZljPL4VnwQLg4ovlgSjUQrG//w58/LH6nU1xobEmFB4GoQQLboYzS0kB3n+fQSgc9eqp7aPz8mT0lmH4cGDvXu47FI3162V5pFA1IECaQn/7jaskhKO6mlD79sDBg/LH6xiEEuzvfw8cYWRWXg7MnAn8+qt6jUGoalWNkjvpJBn2zia5yDVoIEsj+XwScIL7I5OTgR9/lGAFAKefLvvnkFXDhkBWlqRD1YQA1oYABqGE8vtlzkAw85e9pERGJBm4WkLVLr5YLez6ySeq89xokps1i7PSI9WunVx/O3YA48YBTZvK68ZKCGVlMinYcNlliS+jkxjf3z/+CKxZMggpDEIJ5PMBv/wCFBQAGzaoJWb69wcuv1wmAzZuHLgKgLldmQI1aQL06iXpLVuATZtU3vDhQH6+9G9y7JDnAAAgAElEQVRQ5Fq3Bu64Q67Zm24Cxo4FTjhB8swPUgxC1TO+v2VlMmrT0Ly5NCkzCDEIaZGZCXTrBhx7rPw+Zox0pK9dK/MwjA719HSgUSN95XSCqprkjjtObppskovekiWyV9DIkcADD0gT3IoVqjn52GOBY47RWkTb41yhmjEIaWSeC2Tw+4E//5T0UUdxjlBNqltVe/hwee3QocSWyS2mT5fVoE89Vb22caNq9mQtqGYcIVczBiGNtm2TpU6M2f+ADNM2VgBo21ZLsRylQwege3dJf/ttYJPHsGFyLj/7TE/ZnKy0VAbJDBsW+CBkHprNIFSzmkbIMQgxCGm1bZs8KZmXv//jD5VmEAqPuTb06acq3amT7G/DJrnILVokD0RXXaVeKyiQeWwA0KZNYA2JQjN/h4MHJTEICQYhjbZts45+M5riAAahcFW1egIgT/Lz5skN9OefgccekyHc8+YltoxOM22a9Kt166ZemzdPrfJx6aVsKg5HTUEoP1+CvZel6C6Al23bZn2aNAeho45KaHEc65RT5Ml8xw55Ui8oUPMzTj9dhml37SpNdUlJMlTWvKoxBSopkWD+wANsiqut5s1lEnp5eeggBEhfkZcHILEmpNG2bdbJqKwJRc7nU7P2S0uBqVOBF18ETj5ZtncAVF+RMVfDvP0yBZo3DygqCtzj6vBhVXts3FidV6pecrIMdwcCv9sA5woZGIQ0KSiQ7RuCm+PYJxSdK65Q6WeeAe67T4a8V4VBqGrTpknt8uij1WtffCGBCZCN7VLYhhI243u8d6/aCh3gXCEDg5AmxnBN9gnFRq9eat+lvXtlpr+xmkKw5GSe26oUFspowuHDA19nU1z02rRRaXMzMOcKCQYhTULNEQJUEGrcWK0STTVLTlY3x8OHgUcfrToQtW7NJ/mqfPqp9AkNHapeq6gA5syRdHo6cNFFesrmVOYgZHy/S0tllY/69YGlS4FHHgGuuUa2IrnxRj3l1IVfRU1CzREyd5jzST1yV14JTJok6ZUrZXZ/r16yDJJ5WwdzMxMFmjZNFjA1T7JcsUItJdW3r2xNTeEzB6E33pC+trw8NekXAH76Sb7/FRXSROclrAlpEmqO0J49aggsg1DkzjtPapCA7DHUurXcQFu0UDWipCSgSxdtRbS1AweA+fOtTXEffaTSl1+e2DK5gfm7fOCALIVkDkCAfO+NPbGuuy5xZbMDBiFNQs0R4qCE2klNVXOGioqAhQtl9OGKFbLdgM8nT5sclBDaxx9LjXHIEPVaZaWsnAAAderIoASKjLkmdNRRMuijqv7KlBTvBXoGIU22buWghHi48kqVNm6e7doBX38tq24DsmcOWU2bJrXJVq3UaytWqOHtffvy3EUjeGDCe+8FtoAYkpNlQ0GvbZXOIKRJTaslcKJqdC64QN0o58xR+wm1bQt8842srM2OdavcXBmGXV1TnLmGROFr2VLVfHbskNXHn37a+vcqKmRwgtcwCGlgLNXBmlDs1a2rmozy89VaZ4D0Ba1fz40CQ9mwQUa+mZuCgpviBg/WUzanM09YNVZNuPdea7NcRoZs1Og1DEIaVDVHiH1CsRGqSY6q17u3PAQZO6kCUnPcuVPSF13EprjaMJrkcnOldp6SEtgs5/NJTdNrI+MABiEtapoj5PNJRzpF56KL1ByrTz5RIw6paj6ftS+CTXGxE2rCqrlZzu8Hrr468eWyAwYhDULNEQJUEGrRQpo/KDppadLBC8jSSEuWaC2OI1VWym6/gIw6ZFNc7VS1mva998pKH6mpwPnnJ75cdsAgpEGoOULl5TJPCAh8aqLosEmudlauDGyK89qIrVgzf6fNQSglBfjuO+m7rGrYttsxCGnQvTvwt78Fvpabq1Z4NjoxKXr9+6uZ/bNmqYmAFB42xcVWqKV7DO3aycoeXsUgpMENNwCPPx74mvHUCTAIxUJGhgQiQJZIWbZMb3mcJLgpzrxpIEWnus3tvI5ByCYYhGLP/AT/wQf6yuE0q1apzvMLL2RTXCxU1RxHDEK2YQ5CHBkXG4MGqVFyM2aoiatUPTbFxV7LlqoPOLg5zusYhGzCvM8Ia0KxkZGhdlw9cABYsEBveZyATXHxkZKilkNiTSgQg5BNsDkuPq66SqXZJFezFSvUTfLCC4FGjfSWx02MfqE9e2Q/IRIMQjbBIBQfF12kFi6dPVttUU2hmQO1OYBT7Zn7hczfd69jELIJ46KsU0fdNKn2UlNVv0ZJiQQiCq2sDPjwQ0mnpbEpLtaqG6btZQxCNmEEodatZQkVih3znCw2yVVt0SJg3z5JDx4MZGbqLY/bmAccsSakMAjZwJEjwN69kmZTXOz17KmeQhcsUOeaAk2dqtJsios98z5Nxh5NxCBkC+YLkkEo9pKS1E21vFyN/iKluFgWewVkXlC/fnrL40YMQqExCNkAByXEn7lJzvzET+LTT4FDhyR95ZWywC7Flvm7zeY4hUHIBhiE4u/EE2XpfABYvpwdw8Hee0+l2RQXH6wJhcYgZANcLSH+fL7A2tC0afrKYje5ucD8+ZLOzgbOPVdvedwqK0stqssgpDAI2QBrQokxfLhKs0lO+eAD6SsDgGuu8e6WAvHm86naEIOQwiBkA1yyJzE6dwZOP13Sa9cCGzfqLY9dvPOOSl93nb5yeIHx/T5wQOatEYOQLbAmlDjXXqvSb7+trxx2sXEj8MMPks7JAY47Tm953M7cL7R7t75y2AmDkA0YQSgjgxME4+2qq2QVBQB4913VDOVV776r0qwFxR8HJ1gxCNmAEYSys7laQrw1aSKrAQCykKSXV9auqFCj4lJSAvvMKD4YhKwYhDQrKgIKCiTNprjEuOEGlZ4yRVcp9Fu0SD0ADRgANGumtzxewLlCVgxCmnG1hMTr2xdo0ULSc+ao9dK85vXXVXrECH3l8BLWhKwYhDTjyLjES02VociA7OvixeHae/ZIAAZk18+LL9ZbHq9gELJiENKMI+P0MDfJTZ4M+P3aiqLFlClqUMaIEWqwBsUXg5AVg5BmXC1Bj27dgB49JL1+PfDdd3rLk0h+vwRew4036iuL1zRuLHuGAewTMjAIacaakD633KLSkybpK0eiLV0KbNki6d69gU6d9JbHS7hqghWDkGbmIGSuqlP8DR0q63kBsnSNMUrR7cwDEm6+WV85vMr4nu/dK32SXscgpJl51jSDUGJlZKgBCsXF3higsH8/MHOmpJs0AS67TG95vIirJgRiENLMuAgzM4H0dL1l8SJzk9xrr7l/gMK778pOvoCskMB9gxLP3OzOJjkGIe2MINSypd5yeNWJJwYuaurmAQrBAxLYFKcHR8gFYhDS6PBhID9f0sbkSUo8c23olVf0lSPeli9XK4f37Kk2+aPEYhAKxCCk0Z49Ks2akD5XXSVDZwFg+nT33hheflml77hDXzm8zhyEOEybQUgrc6ckg5A+6emqNlRWBkycqLc88fDHH8DHH0u6ZUvgiiv0lsfL2CcUiEFIIwYh+7j9drWj6MSJqvPeLSZMACorJX3bbWrCJCUem+MCMQhpZA5C7BPSq21bVTvIzQWmTdNbnlgqLlZzg1JTgZEj9ZbH65o2VQ88HKLNIKQV+4TsZfRolf6//3PPcO2pU2V+ECB7BvGBR6+kJLVthvke4FUMQhqxOc5ezjxTtrgGgB9/BL7+Wm95YqGyEhg3Tv0+apS+spBiPAjk5rrnYSdaDEIaMQjZi88XWBt66SV9ZYmVTz8FNm2SdM+ewGmn6S0PCSMIlZUBBw/qLYtuDEIamYNQ8+b6ykHK0KHqgeDjj4HNm/WWpzb8fmDsWPX7P/+prywUyNwk6vUmOQYhjYwgZF7enfSqUwe4+25J+/3A88/rLU9tLFsGfPutpLt1ky28yR7MD50MQqSF388le+zqttuABg0k/e67wI4dessTreeeU+mHHpLmRrIH1oQUBiFNioqAkhJJMwjZS1aWWlGgrAx48UW95YnG2rXA/PmSbt8eGDZMa3EoCIOQwiCkCecI2dvo0UBamqQnTnTepMJnnlHp++4DUlL0lYWszN/53Fx95bADBiFNODLO3po3l2Y5QBaaffZZveWJxA8/ADNmSLpFC2DECL3lISvWhBQGIU0YhOzvwQdl4ztA9hravl1vecL1yCOBae5TZT8cmKAwCGnCIGR/zZureUNlZcD//I/e8oRj+XLg888lfdRR3DPIrowVEwAGIQYhTcwXHvuE7Ou++9RIubfeArZs0Vue6vj9wMMPq9+feII7p9pVaqpsrw6wT4hBSBPzhccgZF+NGgH33ivpiorAm7zdLFiglho65hjg2mv1loeqZ3zvWRMiLcxBiKsl2Nvo0ar55KOPZBKo3ZSXA/ffr35/+mmOiLM7IwgVF8uUDa9iENLEHISaNtVXDqpZVlbgkOe77pJakZ1MmKC27s7J4aZ1TsDBCYJBSBMjCDVqxCV7nODvfwdOPlnS69YBb7yhtzxmubnAY4+p38ePl+0CyN44TFvwUtXECEJsinOG5GTZY8gwZgxw4IC+8pg9/DCQny/pESOAM87QWx4KDyesCgYhDUpKgMJCSTMIOUevXrIpHADs3asGLOj0zTfAm29KOisrcNVssjfWhASDkAZ5eSrNIOQsL7wgN3tAhmwb67PpcOgQcN11alO0p57iSEsnYRASDEIacGScc7VpA/z73+r3m28GCgr0lOWBB4DffpP0WWcBd96ppxwUHQ5MEAxCGjAIOdtNNwF9+kh6xw7gH/9IfBkWLZIRcYAsy/P229JvRc7BmpBgENKAQcjZfD7g9dfVunJvvgm8/37iPn/vXhmtZ3j+eaBz58R9PsUGByYIBiENGIScr3174NVX1e8jRwKbNsX/c0tLZQ6QsdFenz5qtW9ylnr1VP8ia0KUUAxC7nD99fIHkEECF18c35uJ3w/cfrtasaFlS2DKFM4JcjIu3cMgpAWDkHu8+irQvbukt24FBg6UgBQPL72kJsnWqwfMng1kZ8fnsygxjO9/fr7sW+VFDEIaMAi5R0YGMG+ejJoDgO+/l620y8tj+zkzZsiK3oY33wROPz22n0GJZ+4XMk/d8BIGIQ2MIJSSAjRsqLcsVHvZ2bKHj7Hlw2efAVdeKZOSY+HddyWwVVbK748+Clx1VWzem/TiCDkXBiG/X7Zivuce4PHHZTKfsaSJXRhBqFkztue7RbduwMcfq3UAZ88G+vWr3dI+fj/w8svS72QEoBEjZJ8gcgdzS4hXR8i57hY4YQKwdCkwbhzw5JPAscdKILILv5/rxrnV+edLLcgYur1smaxovXJl5O+1fz9wzTWyjYSxIsIddwCTJ/PBxU3MO6yyOc4lnn8+MOhcdx0wZ44a0qpbfr5sFQ0wCLlRnz7AkiVqe47ffwd69gRuvVUCS00qKoCpU4Hjj5efhocfBl55hQHIbcxBaO9efeXQyVWX9C+/AH/8IV9gQ3a2tNV/8YW+cplxUIL75eTIAAVj4IDfD7z2GtCuncwn+vJL2cjMUFoqewE9/7xcu1dfDezeLXkNGgAffCD7Gfl8if+3UHyZ9xLzak3IVXsvGutoGRPADJmZEqDsgEHIG9q1A1askL19HnlEhm0XFQGTJsmfpCQJMOnpEnBCbZI3YID8XQ7Ddi/WhFxWEzI6gY02eUP9+tYO4nHjZFhtmzbSfGK2davKGzXK+jmDB6v8YFOmqLxZswLzCguBSy5RvwcHoREj1LH79gXmzZ2r8iZNsn5ux46S16+fNe/++9WxmzcH5q1cqfJCbQOQkyN5OTnWvLFj1bHB/R6bN6s887bThn79JK9jR2vepEnq2LlzA/P27VN5I0ZYj736apVvbJdhmDVL5U2ZYj3WyBs82Jo3apTK37o1MG/JEpU3bpx6PSUFuPtuoG1buSbNTWmVlXJN/vWXNQCdcw6weLGc3zPOkPcdM8Zapt69Ja9rV2ve+PGqTAsXBubt2qXyRo60HjtkiMo/ciQwb+pUlWduLgTk7xp5Q4ZY33fkSJW/a1dg3sKFKm/8eOuxXbtKXu/e1rwxY9SxGzYE5q1Zo/KefNJ67NlnS54x18ssEfeI1avV63l5cs0aeVdfbT02nvcI87WbSK6qCRkLOAYv5FhWZp23UVAgNwDA+kWrqFB5oUY35eWp/GCHDqk8c5MLIM0y5n6B4CC0b5861hgNZSgpUXmh9qP/6y9p1glVuzJudoD1PBw5ovJCrQa9e3fV/9bqzmF5efXnMDdX8kPtKltUpI4NHuZcWanygr+EgDxNGvlGh76huFjlhZpQauS1bWvNM5/D4KBR0zncu1c+r2NH2W5hyRJg7Vrg4EFg5051ndx+u4zq7NRJfv/hB/W+Bw9a33fPHsnPzLTmFRaqY4MnQZqv71D9VNVd3+ZzGHx9AyovVNPS/v1Vn8PDh1Ve8MMDIOepsFANgzc7eFAda/S3GkpLVV6oUbLG9R18/QKJuUeY71V5eXLNGnmhakbxvEcEP7wniquCkFG1Df7POXTIevFmZalmjrp1A/OSk1Veo0ahP6eqJpKMDJWXnh6Y5/PJ5xo3KnNVHACaNFHHBndAp6WpvPr1rZ+bnV31BdaokTo2Jeh/vG5dlRfcjAnI0jDmn2bVncOUlOrPYfPmkh8qCNWvr45NSwvMS0pSeU2aWI9t2lTlB/ehpKervFBfOCMv+P/F+DcY+cEPOTWdw1at5O+0aiVPt+Yn3CeflMVQAVmU1AhAAJCaqt431HyyFi3kxhrqesjMVMfWqxeYZ76+Gze2Hlvd9W0+h8HXN1D9OWzcuOpzWK+eygsVVFu3lhtrqP2SGjZUx6amBubVqaPyQgWwli0lAIY6D4m4R2RlSfkPHpSg4/OpPHN/kSGe94hQ124i+Pz+4OdF59q6VZ42168HTjhBXquslP/0cePssdDjXXfJKCdAmrB69NBbHiLSq0sXYMsWCUZ22TI+kVzVJ9Shgyxpb+732LxZnnRCtSXrYK5im590/H4ZwXfuuYkvEyWWEyZUU2xV9/02ao0HD1qbE73AVc1xgMwuf+cdWTYFkC2YBw8GjjlGb7kMoYLQhx/K+mNFRcD27XrKRYljTKhesEB+HztWAtHs2XrLRfFR0/fb3HS5b1/opm83c1VNCAAefFCqt6NGydDYXbtCj4TSxeiwTUlRbdRDh0oZBw7UVixKILtPqKbYqun77fW5Qq6rCaWmAi++qLsUVTNqQk2bcvKhF9U0odrYn4i8w+tL97iuJmRnfn9gECLvccKEakosr09YZRBKoOJiNWeDQcibIplQTd7g9eY4BqEEMl9goeZRkPtFMqGavMHrzXGu6xOygw8+AGbOtL5uftJdswb45BPg0ksTVy7SL5IJ1eQN5pqQF5vjGITi4KqrQu98OX++rKAMyKx5BiDv6dBBfu7Zo2auV1bKHJFQ6+iR+3m9JsTmuASqaqIqeYcTJlRTYnFgAiWM+SknVBCqrLQ205D7GBOqDXabUE3xUdX3OyNDrU3nxZoQm+MSyPyUY376+fxz2bZ5+XK5CHv2lKXr33gj8WWk+HvwQeChh2RCdYMG9ptQTbFV0/fb55P7wY4d3gxCrlrA1O5GjlT7fPzwA3DyyXrLQ0T2cMopwI8/ykoqpaXemsjO5rgEYp8QEYVi3A/Ky723mC2DUALV1CdERN7k5RFyDEIJZNSEMjKsm7URkXd5ea4Qg1ACcd04IgqFNSGKu8pK2SsEYBAiokAMQhR3Bw6oOQIMQkRkxuY4iruq5ggREbEmRHHH4dlEVBUvL93DIJQgDEJEVBUv7ynEIJQgnCNERFVp3FitksAgRHHBPiEiqkpysgQigM1xFCdsjiOi6hgPp6wJUVwwCBFRdZo0kZ9FRbKIqVcwCCUI+4SIqDpGcxwA7N+vrxyJxiCUIOaakPliIyICVE0IYBCiODCCUOPGsmcIEZGZ+eHUWOLLCxiEEoSLlxJRdcw1IQYhiqnSUqCgQNLmC42IyMDmOIqbAwdUmkGIiEJhcxzFjfmphoMSiCgUNsdR3DAIEVFN2BxHcWN+qmFzHBGFwuY4ihvWhIioJmyOo7hhECKimqSlAXXrSprNcRRTDEJEVBOfT9WGWBOimGKfkP1t3w60a6e7FJGbMgW46CKgRw+gSxfg7ruBPXt0l8p+Xn4ZuOUW3aWomXF/YE2IYoo1IfvKzwdefx3IyQH++EN3aSLzxBPAL78ACxYAq1YBn3wCzJgBdO8O/PST7tLZx+7dwKOPOmNlauP+cPgwUFystyyJwiCUAAxC9jRoEHDJJcDvv+suSeTWrAG+/hp49lm1I+fxxwMTJwK5ucDQoYDfr7eMdvHAA2rFErvz4uAEBqEEMIJQcjKQlaW3LKR8+imwZAkwdiyQkaG7NJGZOBG4/Xbr6wMGyIPOzz8DX32V+HLZzddfA5WVuksRPi/OFWIQSgDjica8jzy53/TpQElJfN77+++B668HPv888PWkJKBzZ0lv2BCfz3aKigrgX/8CHn5Yd0nC58W5QgxCCWA80bApzhv8fmD0aGDWLKBevciOLSoK7++VlcnfnTvXmldREfjTLcI9N4YJEyRQp6fX/HdLSoCTTtIfuNkcRzFXVgYUFkqaQcgb7rxTbmbvvhtZzXfyZODYY4EtW2r+u48+CpxxBnDjjYGvl5UB//2vpE85JfzPtrtIzg0gOxkvWwZceWV4fz8tTT5j4EBg9eroy1lbXtxdlUEoznSuoL1rF9Cxo0yA8/msf5o1A266qfr3eOABaeIJPjYzE3jttcT8O5xkwgTggw+A994D6tQJ/7jJk4GRI6XmdN55Nd9shw2TEXHBgWb6dODQIXn9nHMiLr4tRXpuAOCRR4Cnnorsc3JygOeeAy68UN9gFdaEKObMF1Kia0KtWsmXqagI+Mc/1Os5ORKg8vLkC16d558HDh6UYwCgXz/gxx+ldjdyZPzK7kSbN8t5fukloHXr8I974w2pPU2ZAmzcCJx1FnD++cBvv0X2+ZWVwAsvyACYyZPl4cHpojk3q1YBDRtKzSlSV10lgzv+9jegvDzqYkeNQYhizg7Ds1NT5YtlaNECaNky/OOzsuQGN3QoMG+etJ2T1ciRwNFHA9dcE/4xb74J3H+/nNdrr5Va6/TpwJAhcrON5In8xReB9eulhnryyZGX326iOTeVlcAzz0hNKFrPPw+sWwc8/XT07xEtNsdRzJkvJJ2rJXTooNLbtkV27DvvSMftO+9wdF9V5s0Dli4FxowJvwZSWgpMnSp9F717q9d9PmDcOODee+Wch2PNGuDxx6U5MLifyImiPTeTJgHDh0tzcbTatgVuvRX4978Tv/qEF2tC8FNcTZni98t4Kb//1Vf1laO01O9PSpJyZGSEf1xent/frJnfv2xZ/MpmB+3aybmJ1tln+/2NG/v9hw/HrEhh27PH72/f3u9/443Ef7ad7N3r9196qfX1rVvl//b668N/r19/9ft9Pr//nntiVbrwHDmi7hc9eyb2s3VJ0R0E3U5nn5BZaqr0U+zYIR3XublA8+Y1H3fPPcCllwK9esW/jE61caNMirzlFrUKcqKUlsoIsMcfB264Qb1eUCC18PbtE1senZYulaWXzjsv8HVjaPf8+ZLXuLEMn69O584ysGPiRBmJ2LBhPEpsVacOUL++lNkrNSEGoTizQ5+QoUMHCUKANMnVFIQWLwYWLZLZ91S1Dz6Qn/36Jf6z77hDOtPNAQiQ/7fCQuvrbnb55fIn2JIl0ofUr58McAjXwIES2D79VPqkEqVJE28FIfYJxZld+oSAwKfirVur/7slJdIu/tJLQKNGcS2W482eLT+Dn8Dj7eWX5Yn95ptlJJfxp6wM+OILoGvXxJbHbYyHipkzE/u5xsPq/v3eWP+PQSjO7FQTMgehmgYnPPmkjPQaNiyeJbKHggI1ofjPPyM7dv9+qSm2axddsP76axnt1bq1/OnVC/jhB5W3e3fo4xYvlqbShx6Splbznzp1gP/8x9lBKNrzEkpenvzcuzeyMhx3nKy2sGBBYlfgNh5WKyqcs/BqbTAIxZld+oSAwBFy1dWE1q+XtvAJE+JfJp1uvVWGm7dsqR4WOnUCunUDLr44vPdYuVKeVk84IbLPrqwE7rpLbq6NGkkg27lTtme45BLp4+nVS9Kh3Hln9cvytGiRuH6MWKrteTF7800JJH/7m/z+2Wcy8m3QoPDKkpQEnHiibKuQyCZpr42QY59QnNlpBe1wmuMqK6WD/ZFH4tupfeSI3Lhzc2P7vp07y+Ke4Zg4sfafZyyR07Zt+Mf4/bJSxVtvSZ/O+PEq74IL5LV//lN+79Sp+s91k1icF7O//13+1EaHDvKgsX594ubHBc8V6tgxMZ+rC4NQnJkXL9U9xyac5rhXX5UAcffd1b9XZaUsDnnJJeGvz2VWt65syOZ0xsz9Vq3CP+bJJ+VG27WrzHkJdsUV6mZ7/vm1L6NT2PG8tGkjPxO5SSBrQhRTdlpBu21bqZFVVMh21n5/YGDcsQN47DFg4UIgpZorw++X2lJpqdSYMjKA/v3jX347MvobGjQI7+9v2iT7FwGyxUCo9eXM/VKJHuygi13PS7Nm8jPS/qTa8Np2DgxCcVRWpjoW7RCEUlLkyW77dqnt7NwJZGer/DvvlNrNaadV/z7/+Ic0Lb7+unQQX3qpdOCee258y29HxhyUcLdsuPdeCd4NGkjHeyhLl8rP7GygS5fal9EJ7Hpe0tLkZ35+Yj4P8N7GdgxCcaRzBe2qtG8vQQiQJjkjCM2cKUu/vPde9cf/9ZfUqO69V35v1UomAb7yijeDUFmZ/Kyu5mj47Tc5VwDQt2/VgWvJEvmZyFrQ7t3AqafG7mbbvLnUbsKZvGvn86I7CLEmRLVip+HZhg4d1BPl1q1Az55SW7vrLhnWW79+9cdnZ6sAZGjUSJrxvMjYMO3w4Zr/7uzZat5HVRS6070AAAzvSURBVKPvjhwBvv1W0onsD2rZUh4wdLDzeYnkISNWvLaIKYNQHNkxCIUaIffgg0CPHsDgwVqKlHDRDhAJNXHQWCgznCBkHrVXVa1x5Uq54QLe6Q+y83kxtmev6eEsllgTopgxX0B2ao4zbNsGfPMNMG1aYkf/APYYoh0LxtBsY4BCdYy/06SJTG4Nxailtm2rhiHv2iV7Fbk1KNn5vBiTmJs2jd9nBOPABIoZO9aEzBNWf/lFRrk980z4m7CVl8tSPosWyQ3j0CGZuT94MDBnjvQrhEPnEG1zjWbKFNkyoKBAvvAXXyzDgFu0CO+9jBtiOCstZGTIz+rmfYTq93j7bem0d2sQsvN5MdZaTORcHfPKG2yOo1qxYxAy14S+/lqa4W69Nbxji4tlJNy+fbLPS0aGPI3m5MhIO6d9YZ54Qm5iCxZIE91PP0nH+AcfAF9+CRx/fM3vceKJ8rOmtfgA2XLbWGculG++AVaskHTPnur1hQtlx1S3svN5MYLQ0UfH93PMjIntBQWBg5vcisv2xJEdg1CbNrK+GCCdrZMmhb8J2/XXS9v8zJnq6bVVK1nmpm5d4Oyz41PmeFizRoLws8+qPqLjj5dVFHJzZRfZcBaPPOUUORfr1skE3upce62c+82bVf+G4aefgOuuUwMdjCfvnTtlY7Vwa5hOZOfzsmGD/OzRI76fE8yoDR08mNjP1YFBKI7s2CeUlKT6Me67L/w1zz79FJgxQ24Ywcv5rF0LnHWWGs7qBBMnArffbn19wAB5YPj5Z+Crr2p+n5QU2Xfm0CEZklydDh0k6BUUAE89JUGrtBR44w3ZguDtt4HbbpO/++uvMqn4wQdluRo3s+t52btXmlk7dQq/eTZWjCB04ID7V9JmEIojO9aEAFnUsVu3yIZVv/KK/AxeVXvDBqk5XHBB7MqXCN9/LzW7zz8PfD0pSQY3AOopuCZDh8rP4PcK5b77ZC7W/PkyLProo6VWtmyZNDU9/bSsjP3ccxLs69cHRo4M+5/lWHY8L8ZgiMsui+/nhGIEofJyecBxNd1bu7rZhReqrXoPHtRdmuiVl/v9KSl+f3q6319WFpj30kvy71u1Sk/ZonX88VLu22+35p16quS9+GJ473XwoJybXr1iW0bS64Yb5Dr47rvEf/bll6t7xx9/JP7zE4k1oTiy0wratbFvnzyRnXiiddLel1/KUis5OTJwoaYVF+zi0UeBM84Abrwx8PWyMrVC9SmnhPdeDRrI6s9ff+2ORVlJ5n3NmSPXdU5O4j/fPELO7YMTGITiyOgTssMK2rXRrJlMygzu1/r9d+k3OessCbRffaVmmNvdsGHAqlXWQDN9ujR/nHKK9PWE6777ZHDGq6/Gtpykx/Tp8hD5wAN6Pp9BiGLCTito14bPJ53Aa9eqILNtGzBqlNQCjDlG06dLR7JTVVbKkN/kZGDy5PBHDQIy2OPhh4HXXlNr85EzVVTIdZCTI9tH6OClIMR5QnFitxW0a+uJJ+TfM2CAjGbKzJSmtx9+kMmd118vgxPC3dLAjl58UTYvmzwZOPnkyI9/8EHgo49kbb0ZM2JfPkqMKVOkSXb16sgeRGLJvCuu24dpMwjFifnCcUMQqqqp6YIL5MvqdGvWyNbREyZY+4nCVaeOTLo8/XS5kd1wQyxLSImwYwdw//2yiki4fYLxwJoQ1Zr5wjFfUGQ/ubmyO+z48bHZDnr2bHm/o4+W/jJyhrIy4JprgEGDpFark5eCEPuE4oRByBlKSyVgPP54YAAqKKh6C/SanHWWzHe59VYPzPFwCb9faq6dOwNvvqm7NAxCFAMMQs5wxx3AVVdZm84WLVKLZkaje3fpXzKWNyJ7KykB+vWT/sDkZN2lYRCiGDD3CZk7Gck+Xn5ZnnxvvlnmQRl/ysqAL74AunbVXUJKlPR0WZLKLrwUhNgnFCesCdnb4sWyBUVFhSwHE8qzzya2TEQG84Or24MQa0JxYq4JMQjZz513SgCqSosWrMGSPqmpqimXQ7QpKuanF97M7MdYmofIrho1koEtrAlRVNgcR0S1Yd7Owc0YhOKEAxOIqDaMIHTkiIzecysGoThhTYiIasMrI+QYhOLEuGiSk2UDLiKiSDAIUa0YzXENGzp7Gwci0sMchNw8Qo5BKE6MJxc2xRFRNLwyV4hBKA4qK9WTC4MQEUWDzXEUtcJCWRAR4Mg4IooOgxBFjSPjiKi2GIQoagxCRFRbDEIUNU5UJaLaYhCiqLEmRES1xSHaFDUuXkpEtcUh2hQ1buNARLVVr578ARiEKEJsjiOiWPDCStoMQnHA5jgiigUGIYoKm+OIKBaM+0dxMVBaqrcs8cIgFAdsjiOiWPDCMG0GoTgwXywNGugrBxE5mxeGaTMIxYFxsWRmAikpestCRM7lhWHaDEJxwG0ciCgW2BxHEfP71cXCkXFEVBsMQhSxw4fVKBbWhIioNhiEKGKcI0REscIgRBHjHCEiihUGIYoY5wgRUaxwiDZFjM1xRBQrHKJNEeOGdkQUK+Z7SH6+vnLEE4NQjJkvFAYhIqqNjAwgOVnSbI6jsJgvFC7ZQ0S14fOp+whrQhQW1oSIKJaMIMSaEIXFHIRYEyKi2jIeZvPzZUUWt2EQijEOTCCiWDIeZsvLgZISvWWJBwahGGNNiIhiyfww68YmOQahGDNfJFlZ+spBRO5gfph14+AEBqEYMy6SzEw1tJKIKFrmIMSaENXICEJsiiOiWHD7hFUGoRgznlQ4KIGIYoE1IQpbWZkavcKaEBHFAmtCFDZOVCWiWKtuYMKOHc6fO8QgFENcsoeIYq2qIdqlpUCXLsDnnye+TLHEIBSl4mLgzDOBdevUa5wjRESxVlVNKDUVqKgAtm5NfJliiUEoSo8+CqxaBTz2mHqNqyUQUaxVNTDB5wOaNQNycxNfplhiEIrCb78BL78s6cWL1WAE1oSIKNaqG5jQvDmQl5fY8sQag1CE/H7gjjuAykr5vbgYmDFD0hyYQESxVt3AhGbNGIQ855NPgAULVBBKSgImTJA0ByYQUazVqQOkpUk6eJ4Qm+M8pqhIakE+n3qtslL6hn76ic1xRBQfVW1sx+Y4j3nqKWDPHuu4/JQU4LXXODCBiOLDuJ+EqgkxCHnETz8B48apZjiz8nJgyhRg3z71GmtCRBQrxv2ksFCGZRuaN5f7jvk1p2EQCtM991T/H11YCGzapH5nTYiIYsX8UFtQoNLNmknLjPkB2GlSdBfAKS68UILQ7t1S/Q3VGfjHHyrNmhARxUrwMO1GjSTdrJn8zMuTWpETsSYUpvvukzlBGzdKv1D79sDddwO//gp8842MmjMuguRkID1da3GJyEWqGqZt3HOcPEKOQShKeXlA27ZA586yfM8llwCHD0tew4aBI+iIiGqjqvXjzDUhp2IQikJJCXDokLX6a1wcbIojoliqqiaUlSVryDEIeYzxH248hQDSOWhcHByUQESxVFVNyOeTh2E2x3mM8R9uDkLFxTJUG2BNiIhiy81L9zAIRSFUTYjrxhFRvFS3xTeDkAcZ/+HmPiEu2UNE8VLTStpsjvOY3FwgMxOoW1e9xsVLiShe2BxHAUJNDGNzHBHFS1UDEwDnByGumBCFvLzA/iCANSEiip+aakLG+nGlpXJ/KisDOnVKbBmjxSAUhdxcaxBiTYiI4iUzU4Zj+/0SZMaMAXbtknvRpk3yev36asI8APz8M3DssfrKHC4GoSjk5QEnnBD4GgcmEFG8JCXJxNT8fGDvXuDZZ+U186r+5gDUrJms5uIE7BOKApvjiCjRjPtKaSlw9tlVLw2WnAzccouspOAEDEJRYHMcESWacV/Jz5dNNKtSWQncdFNiyhQLDEIRqmndOIA1ISKKPeO+cuQI0LGj7HGWFHQHT04G+vSRVf6dgkEoQqFWSwBYEyKi+AqesPrYY0CLFoHNchUVwG23Jb5stcEgFKFQ68YBHJhARPEVPEy7fn3g1VdlZJyhWTNg4MDEl602ODouQqGW7AFkg7vBg4GiIud0CBKRc1xxBdC1qwSjJk3ktUsvBfr2BRYulN+dNCDB4PP7zXGUavLOO8D118twSPOyPUREOvz2G9Cli9SItm51Vn8QwJpQxEKtG0dEpEunTsDNNwM//eS8AAQwCEXs0CGgdWvdpSAiUqobsm13bI6LUG6u/OnWTXdJiIicj0GIiIi04RBtIiLShkGIiIi0YRAiIiJtGISIiEgbBiEiItKGQYiIiLRhECIiIm0YhIiISBsGISIi0oZBiIiItGEQIiIibRiEiIhIGwYhIiLShkGIiIi0YRAiIiJtGISIiEgbBiEiItKGQYiIiLRhECIiIm0YhIiISBsGISIi0oZBiIiItGEQIiIibRiEiIhIGwYhIiLShkGIiIi0YRAiIiJtGISIiEgbBiEiItLm/wMPaxGko5HYTAAAAABJRU5ErkJggg==\n",
"text/plain": [
"Graphics object consisting of 17 graphics primitives"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ah = var('ah')\n",
"E_ah = 0.15\n",
"ah_lim=1.08\n",
"ah_plot=plot((1/2)*(ah^2-ah^4), (ah,-ah_lim,ah_lim),aspect_ratio=12)\n",
"ah_plot+=text(r'$E_{\\hat{a}}$', (-1.35,E_ah),fontsize=8)\n",
"ah_plot+=text(r'$0$',(-1.3,0),fontsize=6)\n",
"ah_plot+=text(r'$\\frac{1}{8}$',(-1.32,0.125),fontsize=6,aspect_ratio=1)\n",
"ah_plot+=text(r'$0$',(0,-0.006),fontsize=6,aspect_ratio=1)\n",
"ah_plot+=text(r'$1$',(0.95,-0.006),fontsize=6,aspect_ratio=1)\n",
"ah_plot+=text(r'$-1$',(-1.0,-0.006),fontsize=6,aspect_ratio=1)\n",
"ah_plot+=text(r'$V_{\\hat{a}}=\\frac{1}{2}(\\hat{a}^2-\\hat{a}^4)$',(0,-0.04),fontsize=12)\n",
"ah_plot+= plot(0.125,(ah,-1.25,1.25),linestyle=\":\")\n",
"ah_plot+= plot(0,(ah,-1.25,1.25),linestyle=\":\")\n",
"ah_plot+= plot(E_ah,(ah,-1.25,1.25),linestyle=\":\")\n",
"ah_plot+=arrow((1.11,-0.08),(1.13,-0.08-.02),width=.5,arrowsize=1.5)\n",
"ah_plot+=arrow((-1.13,-0.08-.02),(-1.11,-0.08),width=.5,arrowsize=1.5)\n",
"ah_plot+=arrow((0.97,1/16),(1.00,1/16-.02),width=.5,arrowsize=1.5)\n",
"ah_plot+=arrow((-1.00,1/16-.02),(-0.97,1/16),width=.5,arrowsize=1.5)\n",
"ah_plot+=arrow((-.34,1/16),(-.27,1/16-.017),width=.5,arrowsize=1.5)\n",
"ah_plot+=arrow((.27,1/16-.017),(.34,1/16),width=.5,arrowsize=1.5)\n",
"ah_plot.save('plot_ah.pdf',dpi=200,axes=False)\n",
"show(ah_plot,axes=False,dpi=200,figsize=[3.2,2.4])"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAaEAAAHgCAYAAAAMrNILAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XmYHFW9//HPbNkn+74QE8CEBBJMIksCBkQBWcL1XgiCIKAiIJuXH+CCKCjLFS8oKlzUqLkKAldQMIgS9sWAkLAkJCRAEsi+L5NMMslMpn5/fC2repmZnt5Od9X79TzzdM2cqe4zPd39qbPUqQrP8zwBAOBApesKAADiixACADhDCAEAnCGEAADOEEIAAGcIIQCAM4QQAMAZQggA4AwhBABwhhACADhDCAEAnCGEAADOEEIAAGcIIQCAM4QQAMAZQggA4AwhBABwhhACADhDCAEAnCGEAADOEEIAAGcIIQCAM4QQAMAZQggA4AwhBABwhhACADhDCAEAnKl2XQHkT1OTtGaNtGmT1KWLNHq06xoBcGnVKvtM6NtXGjxY6tTJdY1S0RKKkBUrpOHDpYkTpRtucF0bAK498IB0+OHS/vtLs2a5rk16hFCE9O0bbG/a5K4eAEpD+HMg/PlQSgihCKmtlWpqbJsQAkAIoagqKqR+/WybEAJACKHo/Bfaxo2S57mtCwC3wiHUp4+7erSGEIoYP4T27pV27nRbFwBu+SHUvbvUoYPburSEEIoYvztOoksOiDv/M6BUu+IkQihywi+2jRvd1QOAW/v2SVu22DYhhKJhmjYASdq6NRgXJoRQNIQQAKk8ZsZJhFDkhMeE6I4D4osQghO0hABIhBAcIYQASIQQHKE7DoBECMGR8FnRtISA+CqH1RIkQihyOna0hUwlQgiIM1pCcIZFTAEQQnDGf8Ft2WJnTQOIn3AI9e7trh5tIYQiyA8hzwuW7QAQL34I9eolVVe7rUtrCKEIYpo2gHJYvFQihCKJadpAvDU2Stu32zYhhKKjJQTE2+bNwTYhhKIjhIB4K5eZcRIhFEl0xwHxRgjBKVpCQLwRQnCKEALiLdwDQgih6OiOA+It/L4Pfx6UIkIognr2lCr/+Z+lJQTETziE+vd3V49MEEIRVFkZrJpLCAHxQ0sIzrGIKRBfGzYE24QQnPAHI+vrpd273dYFQHH5LaGOHaVu3dzWpS2EUEQxQw6ILz+E+veXKirc1qUthFBEEUJAPDU3B+/5Uu+KkwihyGKaNhBPW7cG1xEjhOAMLSEgnsppZpxECEUWIQTEUzmdIyQRQpFFdxwQT+XWEirhi76Wruefl667zi4atXChXUa7okKaOFHq2jX4Pc+TGhqk5csTXxhz59rvFhItISCeyukcIYkQysrUqdJLL9n2yJEWMtOnSw88kP73PU+aM0c66yxp5crivDAIISCeyq0lRHdcDjZtkj74wLaPP77l36uokKZMke64w7YHDSp83eiOA+KJMaEYeeEFa+VI0uTJbf/+4Yfbi6KmprD1kqQuXexLIoSAOCm37jhCKAcvvGC3ffpIo0a1/fudO0uDBxe2TmH+UVD4RQkg2uiOi5Hnn7fbI49MvzTG448nfr99uzRkSOHr5fNDaPNmqampeI8LwB0/hDp0kGpr3dYlE4RQlrZvl+bPt+10XXGrVkm33574s27dpM99rvB18/kh5HkWRACir5zWjZMIoay99JKt0STZpINkDzxgs+jCBgyQPv/5wtfNFx6UpEsOiL7m5iCEyqErTiKEsuaPB9XUSB//ePDz5mZp9mzpppukY45xUrV/IYSAeNm2rbzWjZMIoaz540Fdukif+YwFzqRJUo8e0gknSHv22Gy4lvadONFC4uKLC1dHQgiIl3KblCARQlnZtUt6/XXbvuYa6bnn7GvuXPuwP/NM6Ygj7IJS6UydakG0aZM0YULmj3vHHdLQofb13HOJZcuXB2WXX24/C4fQzTcH5clmzgzK/vjHxLIdO4KydF2JF1wQlCePOz32WFD2i1+k7jtypJWdeGJq2TXXBPsuWZJY9vLLQdmtt6buO2mSlU2alFp2663Bvi+/nFi2ZElQds01qfueeKKVjRyZWvaLXwT7PvZYYtnmzUHZBRek7vv5zwflO3Yklv3xj0HZzJmp+/pl06alll1+eVC+fHli2XPPBWV33JG677hxVnbUUallN94Y7DtvXmLZggVB2XXXpe77yU9a2ejRqWU/+1mw7+zZiWVr1wZlF12Uuu8ZZwTle/Yklv3+90HZ73+fWLZnT1B2xhmp93vRRUH52rWJZbNnB2U/+1nqvqNHW9knP5ladt11wb4LFiSWzZsXlN14Y+q+Rx1lZePGpZb99KeJf1tYus+IsGnT0r8WCo0VE7IwZ47U2GjbyV1unTtLRx+dfpWCxsbgHKHXX7cJA5/4ROaPW1cnrV5t28kvsH37grKtW+02HELbtgXlyerrg7JduxLLPC8oS/c3bd4clPtjZL7du4OynTtT9129Wtq7N/0JdVu3Bvsmz+zbsycoq6tL3Xfdupb/1taew6am1OcwbMMGK+/QIbVs585g3+Qr2TY3B2XpJohs2hSU++ed+XbtCsrq61P39cuGDUstCz+HfheNr63ncO1aq1enTqll27cH++7dm1jW2BiUbduWuu/69VaebtbWjh3Bvg0NiWXh1/eWLan7btzY8v88/Bwmv76loCzd+XRbtrT8HDY0BGXJBw+StGaN/bxHj9Sy8PvR/yzx7d0blG3fnrqv//pOfv0m/w3JV1RN9xmRvG+610KhEUJZCHfFHXZYavnixbaMT9jWrdJ//If0zDPBffTvn/6IsCXduwdTvJNbWVVVQVmvXnYb/nAPlyfr2jUo809w9VVUBGXhpYB8ffoE5ZVJ7erOnYOydJcYHjKk5RDq1SvYtzrpVdqxY1DWvXvqvgMHJt6GtfYcVlenPodh/ftbeboQ6tYt2Ldz58SyysqgrE+f1H379g3Kk2czdekSlIXXJfT5Zem6XsLPYVVVYllbz+GgQfY76Z7DHj2CfZOfi5qaoKxnz9R9BwywD9Z0r4fa2mDf5PALv357907dt1+/ll/f4ecw+fUttf4c9u7d8nPYqVNQli5UBw+2g5MBA1LLevYM9k0+eb1Dh6AsXYANHGgBmO55CB8IJr9f031GhPXrl/61UGgVnpd87IW2TJ1qExM+9SnpySdTy3fuTH2T/eAH9k/+4hft+09/2l5gDz1UuHquXRucHDttmvToo4V7LADu3XSTdP31tv3II9Jpp7mtTyZoCbXTnj3Sq6/adkuz35IDaNcu6b77pNdes++bmmw84oILpKuusp8tWmTjCvvtl7+6ho+EmJgARF+5rRsnMTGh3V59NeivznQK9s03W1ec3/0zd671769ZI/3whzYYOGFC+kHrXNTUBE12QgiIvnJbN04ihNqtrfGgZH/6k3TnndJXvxr87IUXrJVy331BP/Pw4dKzz6Yf4MwF68cB8cEU7Yhrbg6m3x55ZOurYS9cKF16qXT66bZUT/gF8fzzNs0yPPi6cqXNjEqeWZUrP4R27kw/MwhAdPghVFPjZpJBNhgTysD559s5JCtWWBeaZN1yRx+dOGOmsdE+7D/4IJjqWFUl/ed/Jt7fyy+nnkPxyis2cyrfRy/hfuGNG63FBSCaym3dOIkQyki6kwSztXGjTdcOL/VTX29r0Z1/fv5fOMmrJhBCQDR5XvmtGyfRHVd0Xbta0ITPv7j3Xmsxff3r+X88lu4B4mHbtuDEbkIILerSRTruOOmdd+z7Zcukb35TmjFDGjEi/49HCAHxUI6TEiS645yYMUO6+mqbDffBB3bCarr1pfKBEALiYf36YDvdKg2lihByYPhw6Q9/KM5jEUJAPJRrCNEdF3GEEBAPhBBKUvIUbQDRFA6hclmyRyKEIq9nz2AValpCQHTREkJJqqhg6R4gDgghlKxwCHHhDiCawgeZdMehpPgvyMbG9FdqBFD+/JZQr17pL7xYqgihGGCGHBB9fgiVU1ecRAjFAiEERFt9vX1JhBBKECEERFu5TkqQCKFYIISAaCOEUNIIISDaCCGUNEIIiDZCCCUtHELhFyuAaCCEUNJoCQHRVq7rxkmEUCx07ix1727b69a5rQuA/AsfXNISQknyLydOCAHRQ3ccSp4fQnV10q5dbusCIL/8EKqttZ6PckIIxYQfQhKTE4CoKdcleyRCKDYIISCaGhqChYkJIZSs8IuTcSEgOsp5UoJECMVGuCVECAHRUc6TEiRCKDYIISCaCCGUBUIIiCZCCGWBEAKiiRBCWejXL9hmdhwQHUxMQFmoqZH69rVtWkJAdNASQtkIL93jeW7rAiA/ynnxUokQihU/hBoabPkeAOXPD6HOnaVu3dzWJRuEUIwwOQGIHv+9PHCgVFHhti7ZIIRihKV7gGjZs0fassW2Bw1yW5dsEUIxwtI9QLSEDybDB5nlhBCKEbrjgGhZuzbYpiWEkkcIAdFCCKGsEEJAtITfx4QQSh4TE4BoCbeEGBNCyevdW6qutm1aQkD5ozsOZaWyMjijmhACyh8hhLLjN9nXr5eam93WBUBu/IPJysrERYrLCSEUM34I7dsnbd6c+/1dc410++253w+A9vNbQv37S1VVbuuSLUIoZvI9OWHRIunFF3O/HwDt09wcvIfLtStOIoRiJ9/TtHv3DpYNAVA8mzdLTU22TQihbOR76R5CCHAjCtOzJUIodmgJAdEQhZlxEiEUO+EQCr+Is+WHEBfJA4orCqslSIRQ7AweHGznK4T27JF27879vgBkju44lKXwEdOaNbnfX+/edkuXHFBcdMehLHXtKvXoYdv5DKGtW3O/LwCZI4RQtvwXLC0hoHyFx4TojkNZ8ceF6uulHTtyu69eveyWEAKKy28J9eghde7sti65IIRiKDw5IdfWUM+edksIAcXlh1A5d8VJhFAs5TOEqqvtSIwQAopn507ryZAIIZShQsyQI4SA4onK9GyJEIqlfLaEJEIIKLaozIyTCKFYIoSA8haV1RIkQiiWCCGgvNEdh7IWPnLK5/pxAIqD7jiUtc6dg6nVtISA8kMIoez5XXJr1uS+AjYhBBTX6tXB9pAh7uqRD4RQTPkhtHu3tH17bvfVu7edt9DYmHu9ALTND6GuXaXu3d3WJVeEUEzlc3ICi5gCxeW/ZwcPlioq3NYlV4RQTBUihOiSAwpv506prs62y70rTiKEYiufM+RYxBQonvB4UPhgslwRQjFFSwgoT+H3Ky0hlK18hhAtIaB4ojQzTiKEYiufIdSpk9Sli7R5c273A6BtdMchEsJLfeTrhFVmxwGFR3ccIqFTp2Ash1UTgPJBdxwiw2/Kr13LqglAuQiHULkv2SMRQrHmh9CePbl3pRFCQHH4PRf9+kkdOritSz4QQjGW72nahBBQWM3NwXs1Cl1xEiEUa/m8zDchBBTexo1SU5NtR2FmnEQIxVr4RRzuZ84Gs+OAwovazDiJEIq18Is4Hy2hrVutuwBAYURtZpxECMXa0KHB9qpVud1Xr142wy7Xy0IAaFnUTlSVCKFYy2cIsX4cUHh0xyFS+veXqqttmxACSh/dcYiUqqqgSU8IAaUvit1x1a4rALeGDZNWrJA2bZIaGmw5n7Z4noXWli3Bl99NcNdd0kMP2SSFjRvtfj/2Meneewv7dwBx4IdQhw5S375u65IvhFDMhceFVq+W9t+/7X1+/Wvpy19OX/b443a54X37gqWAonLEBrgWpct6++iOi7lsJiccfXQwlpRs3z47mS68Ft3ll2dfPwCmoSG4XEqUDuwIoZjLJoQ++lHpqqukygxePYMHSyefnF3dAASiODNOIoRiL9tp2tdfb33SrXUJVFZKl11mEyAA5Cb8/iSEEBnZhlC3btKPf9z6JSAqK6UvfSn7ugEIhN+fw4a5q0e+EUIxl8sJq5/7nDRlSvqWTnW1NH26nYsEIHcrVwbbhBAiY+DAYGwn/CLPREWFTclOt15cU5N06aW51w+ACb8/wweP5Y4Qirnq6uCSDtmcsDp+vHTJJYmtoYoKaexY6cgj81NHALSEEGH+UdX69dLeve3f//vfl2prg+89T7ryyuicxwCUAj+EqqqicVlvHyGEhKZ9Npd06N1buu224PuuXaWzz869XgACfggNHhytGaeEEPKymvYXvyiNHGnb06dbEAHIj927bQksKVpdcRIhBOUnhKqqpB/9yMaYLr44P/UCYKI6PVti7Tgof9cVmjZNamzMvT4AEkU5hGgJIa8XtwOQf1Gdni0RQlDikRUhBJSeqE7PlgghyKZ7+tOpCSGg9BBCiLQOHaQBA2ybEAJKDyGEyPP7mdeutSV3AJQOP4RqaqK3HiMhBElBCDU3WxABKB1+CA0dmtl1vMpJxP4cZIsZckBp2rlT2rbNtqPWFScRQvin4cOD7Q8/dFcPAInCB4VRm54tEUL4p/32C7ZXrHBXDwCJojwpQSKE8E/hEKIlBJQOQgixEO6OoyUElA5CCLEwYICdLyTREgJKCSGEWKisDF7gtISA0kEIITb8caHt2+0LgHt+CHXqJPXp47YuhUAI4V8YFwJKi+cF78Vhw4I1HqOEEMK/MEMOKC1btkj19bb9kY84rUrBEEL4F1pCQGn54INgO/z+jBJCCP9CSwgoLeH3ISGEyKMlBJSWcAjRHYfIC0//pCUEuEd3HGKlU6fg4na0hAD36I5D7PjjQmvWSHv3uq0LEHd+CFVXS4MHu61LoRBCSOAfbXmetHq127oAceeH0NChFkRRRAghATPkgNJQVydt3WrbUe2KkwghJOG6QkBpiMPMOIkQQhKusAqUhjhMSpAIISShJQSUhjhMz5YIISShJQSUBrrjEEu9e0tdutg2IQS4Q3ccYqmiIjjq+vBDqbnZaXWA2PK74yoqonkxOx8hhBQjRtjtnj3SunVu6wLEld8SGjxY6tDBbV0KiRBCCj+EJGn5cnf1iIoPPyy/7pSZM6Xjj5eOOEI68EDpa1+T1q93Xavyk+3zuHu3tGGDbZfba6e9CCGkIITyY/t26Ze/lCZNKq+ZhjfcIL37rvTEE9Irr0iPPCI99JA0bpy0cKHr2pWPXJ7H8OslypMSJEIIaYRDaNkyd/UoZ6eeKp12Wvk9f/PmSS+9JN1yS3Ap6bFjpXvusSPz6dNtSSe0LtfnMS7TsyVCCGmMHBls0xLKzqxZ0nPPSbfeKnXt6ro2mbvnHumrX039+Ukn2czJRYukZ58tfr3KTa7PY1xmxkmEENKgO670PfigjRvk29y50nnnSX/9a+LPKyulAw6w7QUL8v+4UZPr80hLCLHWvbsdrUmEUKnxPOnKK6U//tGu/5RvjY3Szp3SY4+llu3bl3gbB7t3S4ce2v7gzfV5DL/vwj0TUUQIIS2/NbRqFdcVKiWXXWYfiL/7XTDWkE/XXy8dfrj0pS8l/ryxUVq82LYnTMj/45aqzp2lGTOkU06RXn018/1yfR79scSKClpCiCk/hJqbW5/ZtXatHal17GhvmOSvfv2kL3+59ce69lrrpkjet7ZW+vnP8/c3lbu775buv1+6997U80bWrJE+9znpk5+UDjlEGjRI6tGj/a2WM8+0mVzJH5APPijV19vPP/GJ3P6OUrBuXeavrUmTpP/6L+nTn858okmuz6P/OEOH2nsrygghpJXp5IRBg+wNs3On9J//Gfx80iQLqI0b7UiyNbfdJm3bZvtI0oknSm+8Ie3YIV10UfZ/Q5QsWWLP749/3PIVNg8+WBo1Snr/ffuQHTNGqqrK/bGbm6Uf/tDua8YMO2AoZ2+/bQFw662Z73PWWTap4Oyzpaam7B430+exrk7atMm2o94VJxFCaEF7JyfU1Nib1DdggDRwYOaP1727vUmnT5cef9z64RG46CLpox+VzjknffngwdK3vy39z/9IvXrZz447Lj+Pffvt0vz51nL42Mfyc5/F1NRkB0SPPip94Qv22nrvvfbfz223SW+9JX3/+9nVI9PnMdzaIoQQW9nMkAvvE57dk4nf/tYGgX/728KMdZSzxx+Xnn9euu66tlshy5bZB65kXXO5mjdP+u53rSsweXyjHGzfbgdERx0l/ehH9hr96Eezu69hw6SLL5b++7/bv3pEe55HQghQdiG0337Bh2R7QmjTJunqq+0IMer939m49VabrfjZz7b9u/65J506SZMn5/a4GzZIp58u/exn0iWX5HZfrvToIW3eLC1daudt3Xij1L9/9vd36aV2sHTbbZnv097nkRACZDNy/BZJpoOxNTXBeEV9fbD2VVuuukr6t3+Tjj66/fWMurfftjPvTz89s4D2Q2jy5NymcO/da4/53e9KX/xi8PO6uva3cqPkgANsPOmee2wcsy3ZPI+EECD7wBsyxLbbc65Qe7vknnpKevJJ6Qc/aFf1YuP+++32xBMz+30/hHIdD7r0UhuMP//8xJ8/+aS1KOLslFOkXbtsVYy2ZPM8EkLAP/mBsmmTzX7LRHixxbbCa/du62P/8Y+DwXQkevRRuz3mmLZ/d8kSm6ot5TYe9JOf2BH/hRfaoL7/1dgoPf20NHp09vcdBf4BwcMPt/572T6Pfgh17WqnOEQdIYQWZTMuFA6htlpCN95og8RnntnempWPujqbai5JK1e2b98tW2yNseHDMwtpvxVUW2vTs2+4waZsDx5s07fvvrvtxUefesq6R7/xDeteDX916GCz7+IeQmPG2NWHn3ii5RO5s30e9+0L3jcjR8Zjkg4hhBblOkOutX3mz7d+9bvvzq5upe7ii20q8MCBFiaStP/+FgYnn5zZfbz8soXGIYdk9vvPPGO3w4dLU6fa2f5vvGGto8svt66hK65o/T4uu6z1E1wHDJB69sysPlFVWSmNHy81NNhBQjrZPo+rV1tLSYpHV5wkVbuuAEpX+E2Q6eSETLrjmpulr3zFzmspxLVS9uyxD+5MJ0Zk6oADbGHKTNxzT+6P5y/vksmlnT0vGGNYutTGK8LjQhddJN1xh83QOvts6cgjW39MtG7ECDtImD8//Tlt2T6PcRsPkgghtKJQ3XF33WVB8bWvZVuz1nXsaBcTK3dLl9rtoEFt/+7ChbY6hWQBmG5iwvjx9rzMnNlyCCEzQ4fabb4v8hcOof33z+99lypCCC0KH4n5H4htGTbMliXZt8+uieJ5if3aq1ZJ3/mONHu2VM2rr1V+qPTo0fbv+uNBffpI556b/nf8bp7XX8+9bnHnTxjwl9fJlzi2hBgTQosGDbJxBcnWI8tEdXVwlLhnTzBby3fZZXadlY9/PHXfHTtskcgRI2wZn1tuyb7uUeDPSMzkfB9/PGjq1JYHs1etslt/ogSy578vtm/P7/0SQkBIZWXQJbBsWeYrMrfUJffww7Z8yU03pd+vttbOnzj/fPugzMeyM+XMb7m01WJsbrZlfaSWp3Lv3WvjF5INihfLunV2vlm3bvn5GjnSDm5cK3QIxeESDj46RNCqAw+0s/YbG22KcSYTCUaMCD4Uly+XpkyxqcpXXGFTU7t1a33/v//dlqk57LCcq1/WunSx24aG1n9v0SJp61bbnjo1/e+8/nownXjixPzULxMDB9qMr6jJ9AChvfxu7yFDCnPRwlJESwit8i9FLGW+8nC6GXJf/7p0xBHStGmt79vQIL34onTCCeV5yYB011Rq66sltbV221YI+TOxunWzKeDphC8zHV7tHNnxL63e1gFVe8TtEg4+WkJoVTiE3n/fxmzaktwdN2eO9MADmc0keuEF+9D9zGfaW9NAKUzRzgd/arY/QaEl/tHz+PEtB/ef/2y3o0fn7xIPceaPq/Xtm7/7DM9ADc9MjTpCCK1KDqFMhN9A775r5wTdfHPLF2MLmz3bWgfV1bZfRYV1N33729Y6yoTLKdrhFQlmzpR+/3s7wt282U5S/eY3Mx+T8cfj2lppoabGbv21/pLNnSu9+aZt33hjPM7CLzR/kkc+Wyzhnobw+y7qCCG06sADg+1MQyjcEnrpJeuGu/jizPadPduO5t9+2y7tUFEh3XefdeMtWJD9tWCK7YYbbAzmiSfsb1i40EL0/vttJtvYsW3fx/jxdtvWOVrjxtltSwP23/mO3Z52ml00ELnzQyifr8dwCIXfd1FXhr3uKKYhQ4JLCGQ6JjR0aHB0Xl0t/eIXmY3vrFtnQXP22dZy8o/YJ02yD/T//d/219+FefMsfG+5Jfgbxo61k0g3bLAgaGsNN0maMMEWsXzrLZsB15JjjrGuu9deS53BeOedNh40YYL0u99l/SchyYIFdnvEEfm7T0IISCM8TXvp0symaVdWBuMZV1+d+dpns2fbbfLJlv7aa+VyHZt77pG++tXUn590ks36W7QoOLm0NdXVdu2a+nrpnXda/72777YQv+km+x/t3m1db1ddZa3Ip54KJjrE1ebNtiq5HyBr1ki/+Y1diTaTgwLfpk3WRbr//vmd7k4IAS3w+6f37s18uu2YMTZTy+8KysTs2daCmjIl8efz5tltMc9vycXcuXZCbnhGmmTh7D+X/gdhW/zus+T7SnbKKbZe3F/+YicZjx4t/eMf0p/+ZB+8cb5UxiGH2MrVffvaxRP9g5rGRrvQ3ODB9rrr3l1asaLt+/NPP8jkSrft4YfQgAHxOmAghNCmbMaFZs2yD1r/pL5MvPqqdb3558f4nnrKbstlVldjo6128NhjqWV+SzLTE38/+1l7PvzZba056SR7DjdssCWTHn+87SnxcbBggR1AeV7LX01NNoFkv/3avj///5rPS5DU1Unr19t2nFpBEiGEDGRzrlA2VqywsYuwujpbReHAAzO/uqhr118vHX649KUvJf68sTE4pyf572xJjx7Sl79sY0xRWJS13DU02AHBpEn2lS/hgztCCEiSzTTtbHzkI6kn/91xh73xf/5zWxi1HJx5pvTKK6lB8+CDNr4zYYKN9WTq6qttcshdd+W3nmi/Bx+07rxrr83v/YYP7splBmi+EEJoU7FC6MIL7Zo4/iDxrFnSrbdKv/yldOyxhXvcYmhuln74QwvSGTPatxrEsGHSt75lQfzhh4WrI1q3b5/9DydNkv7jP/J733EIuBOpAAAbLklEQVSdlCBxnhAyMGyYDezu3VvYELrySvuQPfVUm0W2Y4d1Q6Vbcbvc3H67LSA6Y4b0sY+1f/+vf136wx+k//f/pIceyn/90LaZM6079dVX87+kVJxDqMLz2jM5EXF10EH2BuzUybqUynFdN1fmzZOOPtqC6JJLsr+f5cttUdcf/tBWGkfxrFplJwV//ev2lW+TJ9uVWiV7fyVPzokyPkqQEb9LrqEh9RpBaNmGDdLpp9tltXMJIMmWQ3r0UeuamzMnP/VD2xobpXPOsRZ6IQJICiadDBkSrwCSCCFkqFjjQlGyd68F0He/a+ej+Orqsj/xdvJk6W9/s2WQ6uvzUk20wvOs1XnAAdKvf12Yx9i61U6kleLXFScRQshQ+M2xZIm7epSTSy+VzjortevsySdtAka2xo2z8aWuXXOpHTKxe7edGjBjRuFmZ8Z5PEgihJChUaOCbUKobT/5iR09X3ihnQjpfzU2Sk8/bSsaoPR16ZK6jFS+xT2EmB2HjIQ/NAmh1j31lK3Ztm+f9I1vpP+dW24pbp1QuuIeQrSEkJHBg4MTSf2z/pHeZZe1vizPgAFSz57Fqw9KW5xPVJUIIWSooiLokvvgg7YvOR1nixe3vk7ZunWua4hS4vcsVFbG67LePkIIGfNDqLmZGXJAPnhe0LMwcqSdhxc3hBAyxrgQkF+rVgVT7eM6WYUQQsbCM+QYFwJyF34fEUJAG2gJAflFCBFCaIcDD7QJChItISAfwu+jgw5yVw+XCCFkrHNnafhw216yJLjkAoDsvPNOsB3u7o4TQgjt4ncZ1NUx1RjIld8S6tdP6tPHbV1cIYTQLizfA+TH9u3S2rW2HdeuOIkQQjuFB08ZFwKyFz6Ii+ukBIkQQjvREgLyIzweRAgBGaIlBOQH07MNIYR2GThQ6t7dtsNHcgDah+nZhhBCu1RUSGPG2PaHH0o7d7qtD1Cu/BDq1Enabz+3dXGJEEK7HXxwsL1okbt6AOWqsTFYBHjUKFtBO65i/KcjW2PHBttvv+2uHkC5ev99u9KuFO/xIIkQQhbCIbRwobt6AOUq/L4Jv5/iiBBCuxFCQG7CPQjh7u04IoTQboMGSb162TbdcUD7EUIBQgjtVlERtIZWr5a2bXNbH6Dc+CHUqVM8L+kdRgghK+EuOWbIAZlraJDee8+2x4yRqqrc1sc1QghZCXchMC4EZG7xYqm52bbj3hUnEULIEtO0geyE3y9xnxknEULIEjPkgOwwKSERIYSs9O8v9e1r24QQkDlCKBEhhKz5b6B166TNm93WBSgX/kFbba00bJjbupQCQghZo0sOaJ8dO6QPPrDtgw+20x3ijhBC1gghoH3CpzPQFWcIIWQt/CaaP99dPYBywXhQKkIIWRs3Lth+6y139QDKBSGUihBC1nr0kEaMsO3584MT8ACkFz5YI4QMIYScjB9vt/X10tKlbusClDLPk95807YHD7bTHEAIIUeHHhps0yUHtGzlSmnrVtsOv2/ijhBCTvyWkBQc5QFIFX5/EEIBQgg5Cb+ZCCGgZYRQeoQQcjJ8uE1QkOiOA1pDCKVHCCEnFRVBl9yqVSzfA7TED6GuXaX993dbl1JCCCFnTE4AWrdtm7R8uW2PHy9V8sn7LzwVyBmTE4DWhVcUoSsuESGEnDE5AWgd40EtI4SQszFjpKoq26Y7DkhFCLWMEELOOnWSDjrIthctkhoa3NYHKDV+CFVWslxPMkIIeTFhgt02NUkLFritC1BK9u4NLnUyapTUubPb+pQaQgh5MWlSsD13rrt6AKVm/nwLIkmaONFtXUoRIYS8IISA9MLvh49/3F09ShUhhLw49NBgcsJrr7mtC1BKwu+H8MEaDCGEvOjcORhwXbhQ2rXLbX2AUuGHUFUVM+PSIYSQN/5RXnMz5wsBkl1ny5+UcPDBUpcubutTiggh5A3jQkCiN98MrjhMV1x6hBDyJjzoyrgQkPg+YFJCeoQQ8ubgg6UOHWyblhBACGWCEELedOwojRtn20uWSHV1busDuOYfjHXsyEoJLSGEkFd+v7fnSW+84bYugEvbtknvvmvb48cHvQRIRAghr8JdDnTJIc7mzQu26YprGSGEvArPAHrlFXf1AFxjpYTMEELIq7FjpW7dbHvOHOuWA+LoH/8ItgmhlhFCyKuqKunww217zRpp5Uq39QFc8Dw7CJOkHj2k0aPd1qeUEULIu8mTg+2XX3ZXD8CV5cul9ett+8gj7TpCSI+nBnkXDiH/aBCIk/DrPvx+QCpCCHnnd8dJtIQQT4RQ5ggh5F2vXtKYMbb9xhusqI348UOoslI67DC3dSl1hBAK4sgj7bapifOFEC91dcEl7sePl2pr3dan1BFCKAgmJyCuXn01WDmbrri2EUIoCCYnIK4YD2ofQggF8dGP2tiQZC0hTlpFXPz978E2IdQ2QggFUVkZvAE3bpQWL3ZbH6AY9u0LlqsaNEgaPtxtfcoBIYSCmTo12H7+eXf1AIpl/vzgEiZTpkgVFW7rUw4IIRRMOISee85ZNYCiCb/OjznGVS3KCyGEgpkwIVjM9PnnGRdC9BFC7UcIoWCqq6WjjrLtdeuk995zWx+gkPbtk154wbb79g1O2EbrCCEUVPhokC45RNn8+XY1Vcle94wHZYYQQkExOQFxET7ICr/u0TpCCAU1caLUtattP/cc40KILsaDskMIoaBqamyqqmQXuVu61G19gEJgPCh7hBAKLtw18cwz7uoBFMpbbwXjQVOnchG79uCpQsF96lPB9pNPuqsHUCjh1/Wxx7qrRzkihFBwEycG68g99ZR1XQBRMnt2sH388e7qUY4IIRRcVVXQGtq2jesLIVrq66WXXrLtj3xEOuAAp9UpO4QQiiJ8dPjEE+7qAeTbCy9Ie/fa9vHHc35QexFCKIpwCIW7LoByR1dcbio8jzM3UBwHHWSXdKiqkjZvlnr0cF0jIHdjx0qLFtmMuE2bgvFPZIaWEIrGP0rct0969lm3dQHyYdUqCyBJOuwwAigbhBCKhnEhRE34dfzpT7urRzkjhFA0U6dKHTrY9mOPsYQPyt+sWcH2SSe5q0c5I4RQNN26SccdZ9urVklvvum2PkAuGhqCk1T797fuOLQfIYSiOvXUYPvPf3ZXDyBXzzwj7dpl2yefzFI92eJpQ1ERQoiKcFdc+HWN9mGKNopu4kTp9ddte+VKaehQt/UB2svzpP32s27lDh3slAP/UvZoH1pCKLpp04Ltxx5zVw8gW2++aQEkSZ/8JAGUC0IIRRcOIbrkUI7Cr9vw6xntR3ccii65K2PTJqm21nWtgMwdfLC0cKFt06WcG1pCKLqKCum002x7797EAV6g1C1aFATQ5MkEUK4IIThxxhnB9gMPuKsH0F5/+EOwPX26u3pEBd1xcGLfPuuSW7NGqqmRNmyQevZ0XSugbXTF5RctIThRVRW0hhobpUcecVsfIBPhrrgpUwigfCCE4MyZZwbbDz7orh5ApsJdceEuZWSP7jg443l2OeQVK6TqamndOqlPH9e1AtLzPLt20Dvv2Pd0xeUHLSE4U1ERDOw2NUkPP+y2PkBr5s4NAujoowmgfCGE4NTnPhdsz5zprBpAm8Kvz/PPd1WL6KE7Dk55njR+vLRggX2/aJFdBhwoJXv2SIMGSVu3Sp07W9dx9+6uaxUNtITgVEWFdMEFwfe0hlCKHnvMAkiS/v3fCaB8IoTg3Dnn2MQESfrf/7Up20ApCR8cnXees2pEEiEE5/r1C67Hsn699Le/ua0PELZmjfTXv9r2kCG2ajbyhxBCSfjiF4PtX/3KXT2AZDNm2Aofkk1IqKpyWp3IYWICSkJTky3js3atXSZ52TJp+HDXtULcNTbauWxr1tjrcvlye50if2gJoSRUV0sXX2zbzc3S3Xe7rQ8g2Qrva9bY9qmnEkCFQEsIJWP9emnYMDv67NXLrjfUpYvrWiHOjjtOeuYZ237iCen4493WJ4poCaFkDBgQnLy6dat0331u64N4e+utIIAOPFD61Kfc1ieqCCGUlMsvD7Z/8hM7mRVw4Qc/CLavuMLGhJB/dMeh5EyeLL38sm3/5S/SSSe5rQ/iZ/ly6YADbHyyb1/pww/pGi4Ush0l59prg+3vf5/WEIrv9tstgCRrnRNAhUNLCCWnudnWk3v7bfv+6ac5QRDFs3GjnR6we7eFz4oVXGKkkGgJoeRUVkrf+lbw/c03u6sL4ue22yyAJOnCCwmgQqMlhJK0b5+tpv3ee/b9M89Ixx7rtk6IvlWrbCZcQ4PUsaP0/vtcN6jQaAmhJFVVSd/+dvD9tdcGffRAoXzvexZAko0FEUCFR0sIJWvfPmnCBGn+fPv+/vsTL4IH5NO770pjxtjrrnt3WzqKrrjCoyWEklVVZf3zvm99yy4uBhTCtdcGC5Vecw0BVCyEEEra8ccHZ6ovXy7993+7rQ+i6fHHpUcfte2BA6Wvfc1tfeKE7jiUvLfekiZOtKPUTp1s6vb++7uuFaKioUEaO9a63yRbLurss93WKU5oCaHkjR8vXXmlbTc0SF/9KiewIn9uvjkIoKlTpbPOclufuKElhLKwc6cNGq9cad//6leJF8KDaW62cbOGBrv1PKlrV/viYmyp/vEPacoUa2VXV0tvvmmtIhQPIYSyMWuWNG2abXftah8YBxzgtk7F1tRk50698459LV5swbx2rX3t2NHyvp06SYMG2bTjYcPs2jhjx0rjxkmjR0sdOhTv7ygFu3ZJH/uYzYqTbHr29de7rVMcEUIoKxdeaJdblqSPf1z6+9+lmhq3dSqkrVttMdc5c+xvffVV+/DMt+pqa2lOniwdfbR9DRuW/8cpJZdcIt1zj20fdpg9v9XVbusUR4QQysrOnXbukL+SwhVXSHfe6bZO+dTcLL3xhs3Wevxx6y7K5B1aW2utnN69rcXTsaPdSlJ9vX3t2GFXCd2yJbO6jB0rnXyydMop0pFHRusDeuZM6YILbLtzZ3vOR41yWqXYIoRQdubOtSP2xkb7/pe/lL78Zbd1yoXn2QzA3/5WeuAB61ZryfDh0uGHW0AcdJB9jRhh3ZOZqq+XVq+Wli61E4H9r3feCc6TSdarl3TiidLpp9ulNfyAK0evvWYtPf+cs1//OggkFB8hhLL0q18FwVNTY9cd+vSn3dapvdassenAv/udtGBB+t8ZO9bOk5oyxVojhVxGZudO6ZVXpBdflGbPbrkV1r279O//btOYjz22vFpI771nz+XGjfb9xRdL//M/busUd4QQytaVV9rVVyXrUpk9WzrqKLd1akt9vfTII9bqeeqp1PXwOnSQTjjBusE+8xmbPODKxo3S3/4mPfaY9MQT0vbtqb/Tv78tpXTuuXYuV0VF8euZqVWrrAX0wQf2/Sc+Ya+Zjh2dViv2CCGUraYmafp06U9/su9ra+0DvtSuPdTcLD33nAXPww9biyPZ5Mn2QT59uo3rlJrGRunZZ239vocfTj8L76CDpC98Qfr850tvUsP771uL8sMP7ftDDpFeeEHq2dNtvUAIoczt2SOddpodqUvWNTdzZmmc8b5okXW13XuvHYUnGzHCgufcc8trqvnu3TZp4v77rZWUvJ5fRYV1033hC9ZtV1vrpp6+efNscsW6dfb9yJHW5Th4sNt6wRBCBeZ50q23Sps22Ztx+XLppz+VevRwXbPo2LVLOvNM+0D0XX21nQlf7HNfNm60D+ff/tY+/JL16GGtnS98wcYmSrn7KhPbtkkPPWR/74svppZ36SJ99rP29x53XPFPmP31r22FDT8oDz7YuuAGDSpuPdAKDwX1s5953vHHB9/fcovnTZvmrj5R1djoeV/5iudZ7NvXxIme9+abhX/sHTs87957Pe+kkzyvujqxDpLnVVV53imneN7//Z/n7d5d+Pq4smyZ533ve553wAGpz4HkeYMGed7VV3ve/PmFr8u6dZ53+umJjz95sudt3lz4x0b7EEIFtt9+9gHlW7XK3hArV7qrU1Q1N3vej37keTU1wQdPZaXnXXSR561end/H2rbN8x580POmT/e8zp3Tf+hOnOh5d97peevX5/exS11zs+e9/LLnXXKJ5/Xqlf65GT/e826/3fPWrs3vY9fVed6NN3pebW3i4116qeft2ZPfx0J+0B1XQO++ayfAvfGGdOihwc979rQTLM87z13douz1121MaMmS4GfV1dIZZ0gXXWQz6NrbLdTUZNOon3nGuv1eesl+lmzoUBuYP/dc1iCTrBvsL3+xsbG//CU4t8tXUSFNmmTnIB1/vG239xwkz7PznH71Kxt/27o1KOvbV7rrLusCRWkihAror3+1E/uWLrXBUN+wYdZHfvPN7uoWdXv2WNB///ups9H69bMgOuIIW6pm2DAbr/ODafNmacMGacUKm1zw9tt2zky6WW2SfdCdcYatvjxlilTJ2vRpbdokPfigBdI//pH+d2pq7IBt4kSbrLH//jZNvVs3m4Yv2VTxNWvsnJ9XX7VZe/7Ctr6qKulLX5Juusn+3yhhbhti0XbffdYVsG5d4s9Hj7auiva6/XbPGzLEvp59NrFs2bKg7LLLUvc99dSgPNlvfhOUPfxwYlldXVB29tmp+55/flC+aVNi2axZQdnPf56674gRVnbCCa391Sbbv33dOs878EDrlkvXLZTt18iRnnfQQZ7Xt2/+//arrw72Xbw4sWzOnKDslltS95040comTkwtu+WWYN85cxLLFi8Oyq6+OnXfE06wshEjUst+/vNg31mzEss2bQrKzj8/8fGuu87zevbM7/+lUyfPO+88zzvnnOBxly1LrNOzzwZlt9+e+vcccoiVTZmSWnbDDcG+c+cmls2fH5R961up+x57rJWNGpVa9tOfBvs+8URi2Zo1QdlXvpK67+mnB+UNDYll990XlN13X2JZQ0NQdvrpqfdbLGV0rnP58Y+sk7t+GhvTd+W0pa7OlluRUqfF7tsXlIW7I3wbNwblyfxlXKTUxTE9LyjbtCl1382bg/LkEy937w7K0rUiVq+W9u61Ex7bku3fPmCAXabZX2vu3/5Nevrp1lebTqdjR1vBe8oUO5l01CibGv7OO1aez79969Zg3+TXyZ49QVldXeq+69a1/H9u7Tlsamr99bNhg5Wnm224c2ew7+7diWXNzUHZ5s3Bz0eNslbKa6/ZbDXJulDfeCN4TjNVUWEz7z77WbuPnj2lc84JHjd5KaK2nsO1a+21nq5bcPv2YN+9exPLGhuDsm3bUvddv97K001Z37Ej2LehIbEs/PpOt+5fa+/tXbtafm9LQZm/goQLhFAB+d0AyR9Q9fXZTdHu3l0aMsS2k8/yrqoKynr1Sl8XvzxZ165BWZcuiWUVFUFZ376p+/bpE5Qnd0N17hyUdeuWuu+QIZmHUL7+9j/9yd7UixfbB+Ajj0hPPmlhe+ih9nu9e1udevaU/uu/rIvo2GNtiZ1i/O29egX7Ji+J07FjUNa9e+q+Awcm3oa19hxWV7f+HPbvb+XpQqhbt2Bfv8vMV1kZlPXpk7pv375B+T332Af0tm0WRA8+KP3mN/b/Gj/ezutpbrb3Tr9+0i9+Yf+bww5LnJ7v/w3+/SYfBLb1HA4aZL+T7jns0SPYN/m5qKkJytKdBDtggIVYutdDbW2wb3L4hV/f6U5kbu293aVLy+9tKShz2WXJmFABLV9uY0Hz59sZ2pK9ibp2le64w5aSB4A4Ywi1gEaMsMHV8CytJUusuV1qS8sAgAuEUIGdd56dTe77zW9sbIFrlwAA3XEF19gofeMb1v/fo4ctoPiTn6Tvd4dbnmfnAX3ve9Lzz7uuTfaiuFRUVP43SEUIAZL+7/9sUc6dO+2ief5y/+XorrukP/85WNT11lvtOkGPPuq2XtmK0v8GqQghIGTmTOmGG8r7g274cOmWW2zlBsmm4Q4daid0FvKieIUWhf8NUjEmBETIu+/aSg/hJYOGDLGuuKefdlcvoCWEEBAhS5fabfL5L7W1FlBAqSGEgAjxVzvo2jXx5926pV8JAXCNEAIiJN9LRQGFxrI9iKT775cefrjt3zvnHFtPLiryvVQUUGiEECLprLPsK25GjLDb9euDdemam209tvDlRIBSQQgBEbFmjfTNb9oioiefbF1wu3ZJf/87S0WhdDEmBIQ0N6d2ZZWTgw+2r9Wr7bIOY8bYslFRWCqq3P83SI+TVQHZVXBnzJBefNGurTJ5sjR6tF0yutw0NtqyUPX1dv2jESPKe6moKP1vkIoQAiJm2TK7LLZkJ6jSDYdSRnccEDHPPmu3nTpZqwEoZYQQEDF+CE2enP4S1UApIYSAiPFD6Ljj3NYDyAQhBETIkiU2VVtiLAjlgRACIsRvBdXW2vTsG26wqdmDB9vU7bvvtgvEAaWCEAIi5Jln7Hb4cGnqVDtx9Y03rHV0+eXSpZdKV1zhto5AGFO0gYjwPGnAADuXpnNnadas1HGhUaPskg5z5khHHummnkAYLSEgIhYutACSpHvuST8xYfx4u505s2jVAlpFCAER4Y8H9ekjnXtu+t9pbLTb118vTp2AthBCQET440FTp0oVFel/Z9Uqu92xozh1AtpCCAER0NwsPf+8bR9zTPrf2btXmj/ftgcMKEq1gDYRQkAELFoUXL576tT0v/P66xZEkjRxYnHqBbSFEAIiYPFiu+3Wzc4HSuevfw22Tzqp8HUCMkEIARGwdKndjh8vVbbwrv7zn+129GiW9EHpIISACKipsdshQ9KXz50rvfmmbd94Y8sTF4BiI4SACBg3zm737Elf/p3v2O1pp0nTpxenTkAmWDEBiICmJmnkSGnfPmnFCqmqKii7807pa1+TJkyQnnvO1pUDSgUtISACqqttcdJ166SbbrIw2r3but6uukqaNk166ikCCKWHlhAQIY8/bitnf/CBrR83dqx08cUWQkApIoQAAM7QHQcAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABnCCEAgDOEEADAGUIIAOAMIQQAcIYQAgA4QwgBAJwhhAAAzhBCAABn/j9QPB5fwYDLmwAAAABJRU5ErkJggg==\n",
"text/plain": [
"Graphics object consisting of 12 graphics primitives"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"b = var('b')\n",
"E_b = 12\n",
"b_plot=plot((1/2)*(b^2-1)^2, (b,-2.5,2.5),aspect_ratio=.5)\n",
"b_plot+=text(r'$E_b$',(-3.1,E_b),fontsize=10)\n",
"b_plot+=text(r'$0$',(-3.0,0),fontsize=6)\n",
"b_plot+=text(r'$V_b=\\frac{1}{2}(b^2-1)^2$',(0,E_b/2),fontsize=12)\n",
"b_plot+=text(r'$b$',(0,-1),fontsize=12)\n",
"b_plot+= plot(0,(b,-2.7,2.7),linestyle=\":\")\n",
"b_plot+= plot(E_b,(b,-2.7,2.7),linestyle=\":\")\n",
"b_plot+=text(r'$1$',(1.0,-0.4),fontsize=6)\n",
"b_plot+=text(r'$0$',(0.0,-0.4),fontsize=6)\n",
"b_plot+=text(r'$-1$',(-1.09,-0.4),fontsize=6)\n",
"b_plot+=arrow((-1.97,E_b-2),(-1.9,E_b-4),width=.5,arrowsize=2)\n",
"b_plot.save('plot_b.pdf',dpi=200,axes=False)\n",
"show(b_plot,axes=False,dpi=200,figsize=[3.2,2.4])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.1 Fundamental constants"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"%display latex\n",
"LE = lambda latex_string: LatexExpr(latex_string);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"declare units as variables"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"s = var('s', domain='positive'); assume(s,'real');\n",
"GeV = var('GeV', domain='positive'); assume(GeV,'real');\n",
"J = var('J', domain='positive'); assume(J,'real');\n",
"m = var('m', domain='positive'); assume(m,'real');\n",
"meters = var('meters', domain='positive'); assume(meters,'real');\n",
"kg = var('kg', domain='positive'); assume(kg,'real');\n",
"K = var('K', domain='positive'); assume(K,'real');\n",
"C = var('C', domain='positive'); assume(C,'real');"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"fundamental constents from NIST 2018"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"scrolled": true
},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"c = 299792458*meters/s"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"e = (1.60217663400000e-19)*C"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hbar = (1.05457181700000e-34)*J*s"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"k_{B} = (1.38064900000000e-23)*J/K"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"G = (6.67430000000000e-11)*m^3/(kg*s^2)"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\kappa = (1.67743454782835e-9)*m^3/(kg*s^2)"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"c = 299792458 * meters * s^(-1);\n",
"e_charge = 1.602176634 * 10^(-19) * C;\n",
"hbar = 1.054571817*10^(-34)*J*s;\n",
"kB = 1.380649*10^(-23)*J*K^(-1);\n",
"G = 6.67430*10^(-11)*m^3*kg^(-1)* s^(-2);\n",
"kappa = N(8*pi)*G;\n",
"#\n",
"pretty_print(LE(r\"c =\"),c);\n",
"pretty_print(LE(r\"e =\"),e_charge);\n",
"pretty_print(LE(r\"\\hbar =\"),hbar);\n",
"pretty_print(LE(r\"k_{B} =\"),kB);\n",
"pretty_print(LE(r\"G =\"),G);\n",
"pretty_print(LE(r\"\\kappa =\"),kappa);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"use c=1 units with unit of energy = GeV"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hbar = (6.58211956547607e-25)*GeV*s"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"k_{B} = (8.61733326214518e-14)*GeV/K"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"G = (4.41583261432942e-63)*s/GeV"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\kappa = (1.10981978405276e-61)*s/GeV"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"m = c^(-1)*meters;\n",
"J = e_charge^(-1) * C * 10^(-9) * GeV\n",
"kg = J*s^2*m^(-2)\n",
"def conv(*args):\n",
" return [arg.subs(m=m,kg=kg,J=J) for arg in args]\n",
"[hbar,kB,G,kappa] = conv(hbar,kB,G,kappa)\n",
"pretty_print(LE(r\"\\hbar =\"),hbar);\n",
"pretty_print(LE(r\"k_{B} =\"),kB);\n",
"pretty_print(LE(r\"G =\"),G);\n",
"pretty_print(LE(r\"\\kappa =\"),kappa);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.2 Standard Model coupling constants from PDG (2018, 2019)"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"G_F = 0.0000116637870000000/GeV^2"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"m_W = 80.3790000000000*GeV"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"m_H = 125.100000000000*GeV"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"GFermi = 1.1663787*10^(-5)*GeV^(-2);\n",
"mW = 80.379*GeV;\n",
"mH = 125.10*GeV;\n",
"#\n",
"pretty_print(LE(r\"G_F =\"),GFermi);\n",
"pretty_print(LE(r\"m_W =\"),mW);\n",
"pretty_print(LE(r\"m_H =\"),mH);"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hbar v = 2^{-1/4} G_F^{-1/2}= 246.219650794137*GeV"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"v ^{-1} = (2.67327142421274e-27)*s"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"g = \\frac{2 m_W}{\\hbar v }= 0.652904833068782"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\lambda = \\frac{m_H}{\\hbar v }= 0.508082923505546"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"hbar_v = N(2^(-1/4))*GFermi^(-1/2)\n",
"pretty_print(LE(r\"\\hbar v = 2^{-1/4} G_F^{-1/2}=\"),hbar_v)\n",
"v = hbar_v/hbar\n",
"pretty_print(LE(r\"v ^{-1} =\"),1/v);\n",
"g = 2*mW/hbar_v;\n",
"pretty_print(LE(r\"g = \\frac{2 m_W}{\\hbar v }=\"),g)\n",
"lambdaH = mH/hbar_v;\n",
"pretty_print(LE(r\"\\lambda = \\frac{m_H}{\\hbar v }=\"),lambdaH)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.3 Gravitational and weak time scales $t_{\\mathrm{grav}}$, $t_{\\mathrm{W}}$"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"t_{\\mathrm{grav}} = (\\hbar\\kappa)^{1/2}=(8\\pi)^{1/2} t_{P}= 5.01325654926200 t_{P}\\\\= (2.70277015574135e-43)*s"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"tgrav = sqrt(kappa*hbar);\n",
"pretty_print(LE(r\"t_{\\mathrm{grav}} = (\\hbar\\kappa)^{1/2}=\\\n",
"(8\\pi)^{1/2} t_{P}=\"),N((8*pi)^(1/2)),LE(r\"t_{P}\\\\=\"),tgrav)"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"t_{\\mathrm{W}}=\\frac{\\hbar}{m_W}= (8.18885475743176e-27)*s"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"tW = hbar/mW;\n",
"pretty_print(LE(r\"t_{\\mathrm{W}}=\\frac{\\hbar}{m_W}=\"),tW)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.4 The scalar field energy density $\\mathcal{E}_{0}$"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\frac{1}{\\hbar}\\mathcal{E}_{0}= 2.84118595562545 \\: t_{\\mathrm{W}}^{-4}"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"E0 = hbar*lambdaH^2*v^4/8;\n",
"ratio1 = tW^4*E0/hbar;\n",
"pretty_print(LE(r\"\\frac{1}{\\hbar}\\mathcal{E}_{0}=\"),\\\n",
" ratio1,LE(r\"\\: t_{\\mathrm{W}}^{-4}\"))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.5 Seesaw time scale $t_{I}$"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"t_{I}= 1.02756853595816 \\:\\frac{t_{\\mathrm{W}}^{2}}{t_{\\mathrm{grav}}} = (2.54945892615748e-10)*s"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"tI = (3/(kappa*E0))^(1/2);\n",
"ratio2 = tI*tgrav/tW^2;\n",
"pretty_print(LE(r\"t_{I}=\"),\n",
" ratio2,LE(r\"\\:\\frac{t_{\\mathrm{W}}^{2}}{t_{\\mathrm{grav}}}\"),\\\n",
" LE(r\"=\"),tI)"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\sqrt{\\frac32}\\frac{g^2}{\\lambda}= 1.02756853595816"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio12 = N(sqrt(3/2))*g^2/lambdaH;\n",
"pretty_print(LE(r\"\\sqrt{\\frac32}\\frac{g^2}{\\lambda}=\"),ratio12);"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"t_{I} c= 0.0764308558042792*meters"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"tI*c;\n",
"pretty_print(LE(r\"t_{I} c=\"),tI*c);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.6 Seesaw ratio $\\epsilon_W$"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\epsilon_{W}= 3.38842679174089e-17"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"epsilonW = (kappa*hbar/(2*g^2*tI^2))^(1/4);\n",
"pretty_print(LE(r\"\\epsilon_{W}=\"),epsilonW)"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\epsilon_{W}= 1.04068084127939 \\:\\left(\\frac{t_{\\mathrm{grav}}}{t_{I}}\\right)^{1/2}"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\epsilon_{W}= 1.02662576744880 \\:\\frac{t_{\\mathrm{grav}}}{t_{\\mathrm{W}}}"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\epsilon_{W}= 1.05492833683428 \\:\\frac{t_{\\mathrm{W}}}{t_{I}}"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio3 = epsilonW*tW/tgrav;\n",
"ratio4 = epsilonW*tI/tW;\n",
"ratio34 = sqrt(ratio3*ratio4);\n",
"pretty_print(LE(r\"\\epsilon_{W}=\"),ratio34,\n",
" LE(r\"\\:\\left(\\frac{t_{\\mathrm{grav}}}{t_{I}}\\right)^{1/2}\"))\n",
"pretty_print(LE(r\"\\epsilon_{W}=\"),ratio3,\n",
" LE(r\"\\:\\frac{t_{\\mathrm{grav}}}{t_{\\mathrm{W}}}\"))\n",
"pretty_print(LE(r\"\\epsilon_{W}=\"),ratio4,\n",
" LE(r\"\\:\\frac{t_{\\mathrm{W}}}{t_{I}}\"))"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\left(\\frac{1}{2g^2}\\right)^{1/4}= 1.04068084127939"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"\n",
"ratio134=(2*g^2)^(-1/4);\n",
"pretty_print(LE(r\"\\left(\\frac{1}{2g^2}\\right)^{1/4}=\"),ratio134)"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\left(\\frac{\\lambda^2}{3g^6}\\right)^{1/4}= 1.02662576744880"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio13=(lambdaH^2/(3*g^6))^(1/4);\n",
"pretty_print(LE(r\"\\left(\\frac{\\lambda^2}{3g^6}\\right)^{1/4}=\"),ratio13)"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\left(\\frac{3g^2}{4\\lambda^2}\\right)^{1/4}= 1.05492833683428"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio14=((3*g^2)/(4*lambdaH^2))^(1/4);\n",
"pretty_print(LE(r\"\\left(\\frac{3g^2}{4\\lambda^2}\\right)^{1/4}=\"),ratio14)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5.7 Units of action for the two oscillators"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\frac{6\\pi^{2}}{g^2}= 138.915667118982"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio5=6*N(pi)^2/g^2;\n",
"pretty_print(LE(r\"\\frac{6\\pi^{2}}{g^2}=\"),ratio5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 7.5 $K(1/\\sqrt{2})$"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"K(1/\\sqrt{2}) = \\frac{\\Gamma(1/4)^{2}}{4 \\pi^{1/2}}= 1.85407467730137"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"K = N(gamma(1/4)^2/(4*pi^(1/2)));\n",
"pretty_print(LE(r\"K(1/\\sqrt{2}) = \\frac{\\Gamma(1/4)^{2}}{4 \\pi^{1/2}}=\"),K)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 7.6 $\\langle \\mathrm{cn}^2 \\rangle$ for $k=1/\\sqrt2$"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\langle \\mathrm{cn}^2 \\rangle=\\frac{2}{\\pi^{2}}\\frac{1}{K^{2}}= 0.456946581044464"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"cn2ave = N(pi/2)/K^2;\n",
"pretty_print(LE(r\"\\langle \\mathrm{cn}^2 \\rangle=\\frac{2}{\\pi^{2}}\\frac{1}{K^{2}}=\"),cn2ave)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 8. Cosmological temperature"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"k_B T = 28.8142120659094*GeV"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"kT = mH/N((6*pi)^(1/2));\n",
"pretty_print(LE(r\"k_B T =\"),kT)"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"T = 3.34375046077032e14*K"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pretty_print(LE(r\"T =\"),kT/kB) "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 9.1 Solution for $\\hat a$ in co-moving time"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\frac{\\epsilon_W}{\\sqrt{2}} = 2.39597956199416e-17"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pretty_print(LE(r\"\\frac{\\epsilon_W}{\\sqrt{2}} =\"),epsilonW/N(sqrt(2))) "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 9.2 $\\hat a_{\\mathrm{EW}}$"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hat a^2_{\\mathrm{EW}} = \\frac{3^{1/2}\\pi}{8 K^2}\\frac{2m_W}{m_H}(2 E_{\\hat a})^{1/2}= 0.254261938075174 \\:(2 E_{\\hat a})^{1/2}"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ratio6 =N(3^(1/2)*pi/(8*K^2))*(2*mW/mH);\n",
"pretty_print(LE(r\"\\hat a^2_{\\mathrm{EW}} =\"),\\\n",
" LE(r\"\\frac{3^{1/2}\\pi}{8 K^2}\\frac{2m_W}{m_H}(2 E_{\\hat a})^{1/2}=\"),\\\n",
" ratio6,LE(r\"\\:(2 E_{\\hat a})^{1/2}\"))"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hat a_{\\mathrm{EW}} = 0.504243927157457 \\:(2 E_{\\hat a})^{1/4}"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pretty_print(LE(r\"\\hat a_{\\mathrm{EW}} =\"),\\\n",
" sqrt(ratio6),\\\n",
" LE(r\"\\:(2 E_{\\hat a})^{1/4}\"))"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"\\hat t_{\\mathrm{EW}} = 0.125799532989201"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"thatEW = asinh(ratio6)/2;\n",
"pretty_print(LE(r\"\\hat t_{\\mathrm{EW}} =\"),thatEW)"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"t_{\\mathrm{EW}} = (3.20720742285761e-11)*s"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"tEW = thatEW * tI;\n",
"pretty_print(LE(r\"t_{\\mathrm{EW}} =\"),tEW)"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
"T_{\\mathrm{EW}} = 0.501068706214232 \\;\\epsilon_a"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"TEW = inverse_jacobi('cn', e^(-thatEW), 0.5);\n",
"pretty_print(LE(r\"T_{\\mathrm{EW}} =\"),TEW,LE(r\"\\;\\epsilon_a\"))"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 8.8",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.15"
}
},
"nbformat": 4,
"nbformat_minor": 1
}