adding sympy

bmorris3 · Jun 6, 2018 · 4239e48 · 4239e48
1 parent 0ab400f
commit 4239e48
Showing 1 changed file with 171 additions and 0 deletions.
diff --git a/radspots/sympy.ipynb b/radspots/sympy.ipynb
@@ -0,0 +1,171 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "IPython console for SymPy 0.7.6.1 (Python 3.5.4-64-bit) (ground types: python)\n"
+     ]
+    }
+   ],
+   "source": [
+    "import sympy\n",
+    "\n",
+    "from sympy import symbols, init_session, sqrt\n",
+    "init_session(quiet=True)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "tau0, p0, b, P, a, Gamma, delta, sigma, sigma_tau_0 = symbols('tau0 p0 b P a Gamma delta sigma, sigma_tau_0', \n",
+    "                                                              positive=True, real=True) "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "T = 2 * tau0 * sqrt(1 - b**2)\n",
+    "tau = 2 * tau0 * p0 / sqrt(1 - b**2)\n",
+    "theta = tau/T\n",
+    "Q = sqrt(Gamma * T) * delta/sigma\n",
+    "sigma_tau = T * sqrt(6 * theta) / Q\n",
+    "sigma_T = T * sqrt(2 * theta) / Q"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$$ p_0 = \\frac{\\tau T}{4 \\tau_0^2}$$ "
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$$ \\sigma_{p_0}^2 = \\left(\\frac{\\partial p_0}{\\partial \\tau}\\right)^2 \\sigma_\\tau^2 + \\left(\\frac{\\partial p_0}{\\partial \\tau}\\frac{\\partial p_0}{\\partial T}\\right) \\mathrm{cov}(\\tau,T)+\\left(\\frac{\\partial p_0}{\\partial T}\\right)^2 \\sigma_T^2 +\\left(\\frac{\\partial p_0}{\\partial \\tau_0}\\right)^2 \\sigma_{\\tau_0}^2 $$\n",
+    "\n",
+    "\n",
+    "$$ \\sigma_{p_0}^2 = \\left(\\frac{T  \\sigma_\\tau}{4  \\tau_0^2} \\right)^2 + \n",
+    "               \\left(\\frac{\\tau  T }{16 \\tau_0^4} b  \\theta^2 T^2\\right) + \n",
+    "               \\left(\\frac{\\tau  \\sigma_T}{4 \\tau_0^2} \\right)^2 + \n",
+    "               \\left(\\frac{\\tau  T }{2 \\tau_0^3}  \\sigma_{\\tau_0}\\right)^2 $$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "sigma_p0_sq = ((T * sigma_tau / 4 / tau0**2)**2 + \n",
+    "               (tau * T / 16 / tau0**4) * b * theta**2 * T**2 + \n",
+    "               (tau * sigma_T / 4 / tau0**2)**2 + \n",
+    "               (tau * T / 2 / tau0**3 * sigma_tau_0)**2)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA8sAAAAwCAYAAAAmTOjoAAAABHNCSVQICAgIfAhkiAAAF6xJREFUeJztnXmcHVWVx79ZgADRBAhgWDNsQYdAkGWEgdiAQUAHFMOooPjYEQYXBBVU6GEVGCSJoIiiPegMiiAogoIBR5ZhFcIue4OGZWCQTQXCMn+cW9PV9Wrf38vv+/m8T3XXcupU1Tl169577rkghBBCCCGEEEIIIYQQQgghhBBCxDGmaQUaZj/ge00rIYQQQgghhBCiP5nZtAI5uQD4cNNKCCGEEEIIIYToP6YDb1Xwq4MngCk1nUsIIYQQQgghRI8wvgQZOwDLAX8rQVadbAA8DzzbtCJCCCGEEEIIIdpFGZXlZem9ijLALOCappUQQoiCTAYWYO/z8cA84LuNaiTEkskw8CLwJvBnYLtGtRFCCNE4k4F9m1YiJ+cBezathOhi+aYVEKLHGIdF94D5z6PASs2pI8QSyzAwsWklhBBCRJK5njG24Al3Aq7w/T8FG2/8eeAW4BXgAWDHguepgm1Rz3Lb2BYlXBMiK28Af3V/L4PNcrCkz3QghBBCCBFkJeDoOk/4hcD/s7HK8sPYWOYNgAuxccHL1qlYAmtgOvYbs4BfAIuw59BpVJtsTAV+TvEGHCGWRCYDd2CV5kMb1kWIJZVHgd9jnQV7NayLEEKIcL4E7FGVcH/X9Tjgs4HtRwKvYxmyPdbFKm6bVqVUDvYChppWogJ2AU4C5mAfzZ1GtcnGz2mXjVTJUdjH1IvAM8ClwEY9IL9qvaum1/VPw6rA9W4peveZ96reaejna1vNLacC9wAzGtQlK1U+l3585v14Tf2GnlE4/XpfDgXuxK7rReAG4AMR+45321cuU4FxwMlY4hiPWcAWgf3+A7g4sG5NRleWPwjcDzwI7F+mkhn4Dr071jotL9M7leX3ATc2rUSNXAHsg72cZmA+8xSwYsvlV6131fS6/mn5NtZgJnr3mfeq3mno52vzcxq9UwZDtc+lH595P15Tv6FnFE6/3pfdgJ2B9bDI5hOBxcDGEft/AasPlsoE4Dbf/4fTPS7uHuCYwLrdsWzZy2E1+QeA1bEkGPfTzMO5F7uZ/UwvVZavBj7dtBI+DsEaeOpqUJmIjTv9px6TnyT3POB/qC5pW9HnVPV9L5PNsGvdL2TbO4C3ub8nAXfTWz1addJLz9xPUb2r9EX5obE8I344EQvHDnYo9BJVPpemnnmUH5RR5veKHcfZcL/TNrsrA9luPM8BB0VsWxX4CzY0t1R+zEhXfXC88gQsBHswsP4a4Afu760Z3fM8F/h4uSomMgUb09vv9EpleX3M+FdL2rEAb6X4Tfbt/z23bpOM55kN/BR4DAuDfxUbw/Zd4sfsT3Xn2ybj+dJSlfw4uZtj06ccXvI5/eR9Th5V3/c45gDfBK7FwoXeAn6UcMzFwJN0Z9vdEliIjVm+k+iCoQk+yYiPRUUSZfVPaK+vVUURvav2RfmhsQ7mg3dgDVbBYWq9RpXPpYlnHucHRW0YesuOo2y4Tk4BrgL+iHWoPQfcDhxLdbM5tM3uykC2G8444GPAa8R3HvyO7k7ewuyJxbqvB+wa2LYFVlm+H8tqPB34IfAE1vMBdlPO9B1zJHBE2UomsDvwk5rP2QS9Uln+Khb+USXeB/dgzG+Cb//bsZd3lnnI98ZeiH/EwjpOwAqDn2CNM+Nijr3AnTNunyJUJT9O7pXA81Sb2C/Pc/JT9X2PYyFmky8B95HuI31Lt1+tWRwLsCZmAy+RrrI8GPPz+2ebfa0qiuhdtS/KD/uTKp9LE888zg+K2jD0lh23wYZfw4bffR/4OlZhusXptQgrP8qmbXZXBrLd0czA6j+vY/c9asyyxxnuvKUyGauFH0x3OMEB7oQfA/6EtfZfgoVce+xBd2U52ENdNXNZMrLF9kpl+WrgVxWfw/sYT8My2Ev8poznuNudY0qEzCi+gTUorZPxfGmpSn6c3A2wysw5JZ/TT97n5FH0vnSw5z2Q8/jtsKiKMU5Gmo90sHfsY7S/sjcGWIDNOnAa6SrLaWmrr1VFEb2r9sWm/RCK+WJb/bBDsfdLUar0lSb8MM4PitowNP9uyWPHTZclEyLWn4jp/62Sz9c2uyuDNthuh/q/heJsd2msQ3dzLNfWs8QnL9vbnXNq3AmzTtPzPPZg1sLivP3MBO7CQrXXwMYof4jRIc+LGF15Xh3rTq+TWWh+5TaxORW06hRgBrAUNt5sJnAR8L+YvV9JtNM97ZYXYKGn6/q2vRpxzBnYMITtgUci9tkac+T5bt9rgRewOcxvxl40UcTJn4K9xJPCX1+lu1BL0ntf7MVXZQRH3ucE6e571fwWS3KYpZII9n5dC0uK55H0DJvgM9j93YfusqIobfS1JCZioV63M9LTHvZbJYfecVTti/LDEdrkh1WVGx5V2HORMimJOD8oYsNJ11QXeew4zIbr5JWI9Re45foh2/rJ7spAttvNa8BDwK1YJPRC4PMxsu51yy3jTpin2/4ywj9IZpLcQ3gz9vBWx17cOwPH59AhL2/HKvJ313jOOpnISOKysZgxzcTGgjzelFIxrIolRHm+aUV8bOaWa2NT8PwKOBfYChsneTXWWhjU+XDsJX8plnnvaKzV73hGxuz7mYdFYQwAf4jR591uORtL5HAZFnq6KfaiuBwb8hB8vknyJwLH+f6fBnwKe+n+0rf+GUYXamn0fh+WLKLKDOd5n1Pa+95WrnfL2VhGSzA//yH2gfA6ZnM/rV+1/+edWFjdPKxhcvuS5bfN15JYBYvI2hAbU36203cPbIjSYifzWSwJTFa946jaF+WH7fTDqsoNqM6e85ZJaYjzg7w2DL1tx2E23Aa8JFN3Btb3m92VgWw32XbHEh9x5t2bdWP2ycXKdA++H4MNzA6OYw5jVywj9kPAgeWqlsjO2Hy+/coA4S1iQ82pFMtMTL/PVXyeLC3757h9n6V73ufzCB8rMQfrHffPbTse+E+3f9DOz8L8ZXvsJe/9whIWfN/JeIHu5Avz3bZ/KyDfY38nKy6HQBq5y2MfinfFyCmDPM8pz32JokN5YZIDpA//nOT2vdm3birmS2DXs4jqMpAnMR5r0b2fkTFag5QXht02X0vDAnfsKYyeQWJN7ONrMd0h5WXYah2+2LQfQnm+OEB7/LBDsWuqstyoy57TlElpSPKDPDYM5dtxWQyQzo7DbLgJjsDKiDOwCIi3sCR5wflv+83uyqANttuh/m+hKNv9OpYzaxrW634yFjWwc4ysVZys0zJr28ecTP1jpEU0/4gZadVTNGX5GL+VaJ3+wW0737dua+wlHZZyfw1GXvxh+gR/gyEyvOQHnwzZtrHbdnkB+R5nun3iQrLSyN3ArbsyRk4ZZH1OkO++RNGhmcoyWDKPuKR4d1BNcpQ0HIe1pG/lWzdIOZXlNvpaErPdcdcRPuzpN4T7XRm2WocvNu2H0ExlGar1ww7FrqmqcqNOe05TJqUhyQ/y2DCUb8dlMUB6O06y4Tp4itH371eMbgyF/rS7MmiD7XZo5lsozHaHsLHMr2LRBQuA9yfIWcads+wx8j3N9STPezhMtCGF/dIWrKIbr7IcN2faMMWfR9qP8aUwJxsmfIjC2k6Of/qzhZjDhr3APSfMWxh5yRsej5C/upP/m5zy/XgtusHW3Kxs5eTEjdEZptgzzfOcipBV36GM8gfI9i5ZhLVYh7EZ6YeZDFPuu25Lp9epgfWD7viileVe9LV/d8d9JGL7hW77jhnlpiHJF4fpLT/Mo/NQBtkD9IYfDiXIq7LcqNOe6yiTesGGs35jDmQ4Ls6G69QZrIL8YSwq6QlGhhJA/9mdxzD572sv2O5QRvkDlG+7SUxw5zwrbqciqcbbxhbYPJvPRmxfFhujc3uCnIfJNi7hiQz7itF4iX/eHrNPnc9jBpZJ71LCnXBtt3zMLTfG5rY7Dwv1CPJ3BXXykjdcESF/mlsWHY8+BruWJ7BxOUX4m1vGJcIo+kyzPqeizKV7nt+ZwG5YIT4c2LawpPNGsSwj99nPipgtHpBSTpm+NR4br/kA8LUMMtPSVl87BJvVYSpwDzak5Frf9m2dvF9HnHcNt3wos8bJJPlir/khtMsXy/LDsq+pynKjLntOWyYl+R/E+0ETNtymb8woGw5Sh85PYxW727By5DxGElT1m915FLmvev+ms90kvPpH2YlIW8fSWHj1G8SH8w7QriQGwloS38Imoa+StD1X3liVqFD9M9z2Xdz/n3D/RyWp+7LbPj+1pqM50B3/lYjtJ7ntH80p32MdJydriGkYqzlZ15UgK4qsz6kKOu4cAyXIGiB9a+pY7KPh4cD6ZbBkWmFhl3UwmfStzXMDx6bxzzb62kexsPADsKRm38Sm7FvLbR+LlUtRMz6s6rZXlYW0al9sgx9CM2HYVfthh/zXVFW5Uac9pymTkvzPI84P2mLDZTJAOjuOsuE2cDt2DVPoT7srg7bYbof6v4XKtN313Dljcyd5ITpZutXb8POYhCWymO2uZVbMtc6iu+WnbJq+L73y83gaS/8fbKlqCi+z4Ioh26YCB2Ep7r0xKMu5ZZj+k7DWRwjP0JsGLwwpTP6K2EfRn+hOWncIFmXxCpbRcduE82zolklhg2nkPom1yE5PkFWErM+pn5iOtX77W2zHYOFOV2O9u03wKpaFM+znRfNc5/6/IYf8Nvra4dh9/y6WdOwwzP4/7bZ777q3ER4O+0W3/js59G2DL8oP2+eHUF25UaU9B0lTJiX5n0ecH8iGq4+EysNqbvkG/Wl3ZSDbLcd2V3DLKqK7WsVYLIPby9iLPooFxFemRTNcTfU9/sFKehQ3u/0eYXQW04nY+K43gJ1867dx+z/J6BfWBOASt+17ubWGW5yMR0P0+bXbtnvgmLStnn4+5mR9NWafLHK98UPrhWwrg6zPqQo6ToeBEmQNOFlperT2cfv+i2/dNlgr60Lfb0YJepXFIKZzkTHLbfO1pbHQtz0Ccs7CpjfxuM0du1dgvzmYnd5H9rk72+KLbfBDaKZnuWo/7JD/mqosN6qw5zCSyqS0/ucR5QdtseEyGSCdHYfZcF1siGVfDjIWOBHT63rf+n6zuzJoi+12qP9bqEzb7ThZRceo9wyXYxccln1yPDZRd9xcW6I8BrBQtLNJdqCvMnpOvCpI8zE+Hhv/cBs259yDWCr5+Vgr/JvAZ0KO+62T/RhwOhZiOuzWXYi9XPPg6bMQ++h5GAur+RaW2OAtLPQ0yE1Yi6efB7GhClG828l7FvgGNkl9Ebkfd/IOjTlnXvI+p7LpUKyA+BDWOj3EyAfsw751UVMUnY8V1k1lu87DIMUry9AuX/NC7IINsMdgCWo8PoTZ5GtYb+NJjGRufYCRsdZ+fkZ0VE6HdvhiW/wQivliW/2wQ75rqrrcyGPPeUgqk9L6n0eYH7TJhouSx46bLEs+hzXOXIVNf3QyFiX6MCONou/y7d9PdlcGbbLdDvV/C5Vpu/PontM7F0VDbevCG7P2iZBt78Eqb/3EmsB/AfdiDzrY0tUk78XS//8AGwMShzdeYI2E/YqQxiY3cfucg2UL/Rk2D92LWM/3QMRxb8c+3IexF/lzWNjLnII6e/p8G1gfu58vYmHrCwjP+pi11dPPV7AC6k3seorIXRrLSnxTwjnzkPc5lU2HYgXEIPHvz+GQYyZhBeQlOc/ZFIOUU1luk695H03BUNVjsY8XP7sB/40lEPkrVpE5mug5LVfCel3WxnoH5jAyD2ZbfLEtfgjFfHGQdvphh3zXVEe5kdWe8xJXJmXxPwj3gzbZcFEGyWbHTZclG2H2tRCrmL6OzQl+C3YtYaHF/WJ3ZdAm2+1Q77dQ2bZ7LeFzUfctXpr2c0K2HQmcUK86lTMVy0IH9hG1iNGhGE3ijStZFcuQl8QvsXEfYgQvzOTADMdkbfWsUu5R7phNC5xXjOYwwgtqUYw8vlakYSotmzq9/HOOyhebp81+2KZyo0ry+J/8YIQ223Cbkd01T5m2OxUbbjI1acewwfJB5tCeSlgct2KtTWHjkmdRb8/yCljyqnULyDgRCzuM4klGBrc/hbXOhbXGVc2FWMIDP950FX8mXej70TSXxbeteElakqY6CyPYSzcmZF0essg9A5ua5LgSzitsmoSjgIuoPlHhkkYeX3sNS4I0O7B+Ntb7UQYbYwling7ZJl9shrb7YRvLjSrI43/yA6PtNtxmZHfNUrbt7o116EVlWs9EbDrtlrGA7pb4MdiY2LJDNeI4jfxZWT0uw5wsDZuRnMk4D7OAXzAy1qkTss8MLBxykm/dblgP/09In1TtR9hc2cK4HmvBzJK0oqrerrxyZ2HhSb3Q2NZ23omFK01rVo2+JI+vgSVFeg0LL38nNvbpZUbmtyzK6ViZ5ke+2Cxt98M2lRtVk8f/5Aftt+G2I7trjjJtdyngRmzYU2EmEV5BaivHYJU6/0t/EyxrXF0sh/WoblNQziKs1SOJFbFJ0bfOIHsIM7gkdsESKczBeu07EfvdSvEEBitig/bHFZTTD4zFXr55GkBuonsowgPEJ/hqUq4QTVLE18Cm2xnGps76PeXOuHARVmEJIl8UYbSx3KiaKv1PiChkd73P14APlyXs48AqZQmrgfdileUzfesOIzqjZRXMwTJvjykgY2XsOvbCsgX+BQu53jKw3zJYeHnWEOYh0lWW/bxMdGX5GMqZeH064ZmYRXqq6u2quhdNCDGai7HInjUZPWRKvijKRjYlhFhSmIbNz10anw38P0R4trILQ45dBUvx/iTW+vIU2StoWZkAvMLoNOAXYGHBdTGP4pOA74jd19uwBoDpWBKsx7CU8WCV8fPJd0+HchwXV1neCStol3X/F8mQXsZceUs6VbV6qjVViPrYAXvnL6Y794N8UZSNbEoIsSSwbPIu6Vke2Dewbgib2+wdgd/kkOMvwsZbzcBSm28CvL9MBSO4Bksw5SW7WkS9ia8uIV0G6Di+iFU+p/nWvQurcE53/2+DXedC329GSvlDlFtZ3tjp5iU02xmbs/ES7INPCCGEEEIIIXqWdYG/9/0/B5tSwM8Q1sOZhkuxibIPxubS9fMzonsfOxl0DuMEJ2c3YAPgroLysnIF8J2Q9Z5ecb8Bt+/5wI8Dx09z+7yL7ByNVXa932KsMu5fl5SGPa6yvL7TbaPA+hUIvxdCCCGEEEII0Wq8kN6lGKlcHuzWrQY8kVPuOKwH+l+xaS8OwuZB/mdsPO8BWMjPMsAj2HgZb8zrcznP6fE7bFLxWdjY3zqnjAKbwmmFkPVzsYzPcTzuljPpzqa9OVZhfSiHTmdj4egep2A97vN96xblkOvh9dw/E1h/lDu3EEIIIYQQQvQ0N2AV3QlYhTbIEOl6lucBewbW3U13ZsVN6Z7qyeODwP1Y7/T+Kc7psTzWc3orcB5WEa+TIyg2jdNy2DQOR/nWjcEaE75dQK6fIcoNw96P7sr28XTPRSeEEEIIIYQQPcmpwHZYL/NaIduHSK4sz8AqVsGEJNfRPZb3U9gcyEHGY9MWrI7Nj3w/2cYd34hVOJ8GpmY4rgxmAG+Qf+6u92CV/fvc3+thPdKPA1PKUJD0leWJWC/3TGzqqGPc30HbGALO9f2/LzaG+mzgwEKaCiGEEEIIIUQL2AKb1/GwiO1DJFeWjwVuD6xbCnge+FJg/elYErAgW2PTZnjMJduUQqdiPdYPZjimTG4g/7zDBwN/AHbHKsh/xZKllVnpHyJdZXmA8LHVQ759JgAvYBV7IYQQQgghhOhbbiS6N3CI5MryfKxX2M8BWE9vsEfyIqxyHmQOo+dKPhILb07LB7BK3blJO1bETlhv+LiGzl8nh1J8qiwhhBBCCCGEaBVjQ9ZdC1xVQOaNWHbkIxlJGDYX61V+PLDvWGBN9/PrMiZEbtp5esFCvt+k/uReHr/GGgHWaOj8dbKY6EgEIYQQQgghhOgblorZNkRyz/IY4MtYluuXsMrzRyL23QF4DKtw+cc4h4VhBxOGJXEw4YnDhBBCCCGEEEKInmQ8Nt7Yn+Arb8IsIYQQQgghhBCib9gVG/v8EMqoLIQQQgghhBBCCCGEEEIIIYQQQgghhBBCCCGEEEIIIYQQQgghhBBCCCGEEEIIIXqE/wOoTwgdBr9siQAAAABJRU5ErkJggg==\n",
+      "text/latex": [
+       "$$\\frac{\\sqrt{p_{0}}}{\\sqrt{\\Gamma} \\delta \\tau_{0}} \\sqrt{- \\frac{1}{\\left(- b^{2} + 1\\right)^{\\frac{5}{2}}} \\left(- \\Gamma b \\delta^{2} p_{0}^{2} \\tau_{0}^{2} \\left(- b^{2} + 1\\right)^{\\frac{3}{2}} - 4 \\Gamma \\delta^{2} p_{0} \\sigma_{\\tau 0}^{2} \\left(- b^{2} + 1\\right)^{\\frac{5}{2}} + p_{0}^{2} \\sigma^{2} \\tau_{0} \\left(b^{2} - 1\\right) - 3 \\sigma^{2} \\tau_{0} \\left(- b^{2} + 1\\right)^{3}\\right)}$$"
+      ],
+      "text/plain": [
+       "               _______________________________________________________________\n",
+       "              ╱  ⎛                           3/2                              \n",
+       "             ╱   ⎜       2   2   2 ⎛   2    ⎞           2         2 ⎛   2    ⎞\n",
+       "  ____      ╱   -⎝- Γ⋅b⋅δ ⋅p₀ ⋅τ₀ ⋅⎝- b  + 1⎠    - 4⋅Γ⋅δ ⋅p₀⋅σ_τ_0 ⋅⎝- b  + 1⎠\n",
+       "╲╱ p₀ ⋅    ╱    ──────────────────────────────────────────────────────────────\n",
+       "          ╱                                                               5/2 \n",
+       "         ╱                                                      ⎛   2    ⎞    \n",
+       "       ╲╱                                                       ⎝- b  + 1⎠    \n",
+       "──────────────────────────────────────────────────────────────────────────────\n",
+       "                                                             ___              \n",
+       "                                                           ╲╱ Γ ⋅δ⋅τ₀         \n",
+       "\n",
+       "_________________________________________________\n",
+       "5/2                                          3⎞  \n",
+       "        2  2    ⎛ 2    ⎞      2    ⎛   2    ⎞ ⎟  \n",
+       "    + p₀ ⋅σ ⋅τ₀⋅⎝b  - 1⎠ - 3⋅σ ⋅τ₀⋅⎝- b  + 1⎠ ⎠  \n",
+       "──────────────────────────────────────────────── \n",
+       "                                                 \n",
+       "                                                 \n",
+       "                                                 \n",
+       "─────────────────────────────────────────────────\n",
+       "                                                 \n",
+       "                                                 "
+      ]
+     },
+     "execution_count": 5,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "sqrt(sigma_p0_sq).simplify()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "sqrt(p0)*sqrt(-(-Gamma*b*delta**2*p0**2*tau0**2*(-b**2 + 1)**(3/2) - 4*Gamma*delta**2*p0*sigma_tau_0**2*(-b**2 + 1)**(5/2) + p0**2*sigma**2*tau0*(b**2 - 1) - 3*sigma**2*tau0*(-b**2 + 1)**3)/(-b**2 + 1)**(5/2))/(sqrt(Gamma)*delta*tau0)\n"
+     ]
+    }
+   ],
+   "source": [
+    "print(sqrt(sigma_p0_sq).simplify())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.5.4"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}