{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\def\\CC{\\bf C}\n",
    "\\def\\QQ{\\bf Q}\n",
    "\\def\\RR{\\bf R}\n",
    "\\def\\ZZ{\\bf Z}\n",
    "\\def\\NN{\\bf N}\n",
    "$$\n",
    "# Sage Introduction (Sage days 100, Bonn 2019)\n",
    "\n",
    "[Sage (or SageMath)](http://sagemath.org) is an open source software\n",
    "for mathematics which interfaces many softwares and\n",
    "libraries. We are demonstrating its usage in a Jupyter notebook.\n",
    "\n",
    "This worksheet is available in ipynb and pdf format on\n",
    "[https://wiki.sagemath.org/days100](https://wiki.sagemath.org/days100)\n",
    "\n",
    "### Python is an expressive langage\n",
    "\n",
    "$\\Big\\{17n\\ \\Big|\\ n \\in \\{0,1,\\ldots, 9\\}\\text{ and }n\\text{ is odd}\\Big\\}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "S = {17*n for n in range(10) if n%2 == 1}\n",
    "S"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "124 in S"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sum(S)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "{3*i for i in S}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sage adds some mathematical objects and functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "8324074213.factor()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = matrix(ZZ, 3, 3, [0,3,-2,1,4,3,0,0,1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m.eigenvalues()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m.inverse()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As in mathematics, the base ring on which an object is defined matters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R.<x> = PolynomialRing(ZZ, 'x')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P = 6*x^4 + 6*x^3 - 6*x^2 - 12*x - 12\n",
    "P.factor()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P2 = P.change_ring(QQ)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P2.factor()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P3 = P.change_ring(AA)     # AA = field of real algebraic numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P3.factor()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P4 = P.change_ring(QQbar)  # QQbar = field of complex algebraic numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P4.factor()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In sage this concept of base ring or more generally ground set is called\n",
    "a \"parent\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R.<x> = PolynomialRing(ZZ, 'x')\n",
    "p = x^4 + 1\n",
    "p.parent()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Partition([3,2,1]).parent()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Partition([3,2,1]).parent() == Partitions()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Partitions(10^6).cardinality()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Object oriented, autocompletion, documentation, sources\n",
    "\n",
    "Python is an object-oriented language. That means that functions are\n",
    "attached to objects on which they act. For example you would write\n",
    "`14.factor()` rather than `factor(14)`.\n",
    "\n",
    "The SageMath documentation is useful but it is much more efficient to\n",
    "use help directly from the notebook. The three following features are\n",
    "essential.\n",
    "\n",
    "-   autocomplation with the `<TAB>` key\n",
    "-   acces to the documentation with `?`\n",
    "-   acces to the source code with `??`\n",
    "\n",
    "Computing the integral of a symbolic function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(x) = sin(x)^2 -sin(x)\n",
    "f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f.in"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculator\n",
    "\n",
    "Integration (symbolic, numeric and certified numeric)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integral(e^(-x^2), x, -Infinity, Infinity)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integral(1/sqrt(1+x^3), x, 0, 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "numerical_integral(1/sqrt(1+x^3), 0, 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R = ComplexBallField(128)\n",
    "R.integral(lambda x,_: 1/(1+x^3).sqrt(), 0, 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R = ComplexBallField(1024)\n",
    "R.integral(lambda x,_: 1/(1+x^3).sqrt(), 0, 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Roots:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(x) = x^5 - 1/3*x^2 - 7*sin(2*x) + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot(f, xmin=-2, xmax=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r1 = find_root(f,-2,-1)\n",
    "r1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r2 = find_root(f,0,1)\n",
    "r2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r3 = find_root(f,1,2)\n",
    "r3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot(f, xmin=-2, xmax=2) + point2d([(r1,0),(r2,0),(r3,0)], pointsize=50, color='red')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Latex:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M = Matrix(QQ, [[1,2,3],[4,5,6],[7,8,9]]); M"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "latex(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M.parent()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "latex(M.parent())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Graphics:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x, y = SR.var('x,y')\n",
    "plot3d(sin(x-y)*y*cos(x), (x,-3,3), (y,-3,3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Interaction:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('x')\n",
    "@interact\n",
    "def g(f=sin(x)-cos(x)^2, c=0.0, n=(1..30),\n",
    "        xinterval=range_slider(-10, 10, 1, default=(-8,8), label=\"x-interval\"),\n",
    "        yinterval=range_slider(-50, 50, 1, default=(-3,3), label=\"y-interval\")):\n",
    "    x0 = c\n",
    "    degree = n\n",
    "    xmin,xmax = xinterval\n",
    "    ymin,ymax = yinterval\n",
    "    p   = plot(f, xmin, xmax, thickness=4)\n",
    "    dot = point((x0,f(x=x0)),pointsize=80,rgbcolor=(1,0,0))\n",
    "    ft = f.taylor(x,x0,degree)\n",
    "    pt = plot(ft, xmin, xmax, color='red', thickness=2, fill=f)\n",
    "    show(dot + p + pt, ymin=ymin, ymax=ymax, xmin=xmin, xmax=xmax)\n",
    "    pretty_print(html('$f(x)\\\\;=\\\\;%s$'%latex(f)))\n",
    "    pretty_print(html('$P_{%s}(x)\\\\;=\\\\;%s+R_{%s}(x)$'%(degree,latex(ft),degree)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Distribution of polynomial roots"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "m = random_matrix(RDF,700)\n",
    "e = m.eigenvalues()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "point2d([(i, v.abs()) for i,v in enumerate(e)], pointsize=3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from itertools import product\n",
    "R = PolynomialRing(CDF, 'x')\n",
    "roots = []\n",
    "for p in product([-1,1], repeat=13):\n",
    "    p = R(list(p + (1,)))\n",
    "    roots.extend(p.roots(multiplicities=False))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "point2d(roots, pointsize=1, aspect_ratio=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Modular forms\n",
    "\n",
    "For more on modular forms see "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "G = Gamma0(10)\n",
    "G.farey_symbol().fundamental_domain()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M = ModularForms(G, 2)\n",
    "M"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M.basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Polytopes over algebraic numbers\n",
    "\n",
    "(This needs the optional package PyNormaliz.) For more on polytopes, look at the [thematic tutorial on polytopes](http://doc.sagemath.org/html/en/thematic_tutorials/geometry/polyhedra_tutorial.html)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P = polytopes.six_hundred_cell(exact=True, backend='normaliz')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P.plot(point={'color':'red'}, line={'color':'blue'}, viewer='threejs')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Differential geometry: plotting the Ads space $\\RR^{2,3}$\n",
    "\n",
    "For more on differential geometry see [https://sagemanifolds.obspm.fr/](https://sagemanifolds.obspm.fr/)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M = Manifold(4, 'M', r'\\mathcal{M}', structure='Lorentzian')\n",
    "M0 = M.open_subset('M_0', r'\\mathcal{M}_0' )\n",
    "X_hyp.<ta,rh,th,ph> = M0.chart(r'ta:\\tau rh:(0,+oo):\\rho th:(0,pi):\\theta ph:(0,2*pi):\\phi')\n",
    "R23 = Manifold(5, 'R23', r'\\mathbb{R}^{2,3}', structure='pseudo-Riemannian', signature=1,  metric_name='h')\n",
    "X23.<U,V,X,Y,Z> = R23.chart()\n",
    "h = R23.metric()\n",
    "h[0,0], h[1,1], h[2,2], h[3,3], h[4,4] = -1, -1, 1, 1, 1\n",
    "# writing coordinates\n",
    "l = SR.var('l', latex_name=r'\\ell', domain='real')\n",
    "assume(l>0)\n",
    "Phi = M.diff_map(R23, [l*cosh(rh)*cos(ta/l),\n",
    "                      l*cosh(rh)*sin(ta/l),\n",
    "                      l*sinh(rh)*sin(th)*cos(ph),\n",
    "                      l*sinh(rh)*sin(th)*sin(ph),\n",
    "                      l*sinh(rh)*cos(th)],\n",
    "                 name='Phi', latex_name=r'\\Phi')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "graph_hyp = X_hyp.plot(X23, mapping=Phi, ambient_coords=(V,X,U), fixed_coords={th:pi/2, ph:0}, \n",
    "                    ranges={ta:(0,2*pi), rh:(0,2)}, number_values=9, \n",
    "                    color={ta:'red', rh:'grey'}, thickness=2, parameters={l:1}, \n",
    "                    label_axes=False)  # phi = 0 => X > 0 part\n",
    "graph_hyp += X_hyp.plot(X23, mapping=Phi, ambient_coords=(V,X,U), fixed_coords={th:pi/2, ph:pi},\n",
    "                    ranges={ta:(0,2*pi), rh:(0,2)}, number_values=9, \n",
    "                    color={ta:'red', rh:'grey'}, thickness=2, parameters={l:1}, \n",
    "                    label_axes=False)  # phi = pi => X < 0 part\n",
    "show(graph_hyp, aspect_ratio=1, viewer='threejs', online=True,\n",
    "     axes_labels=['V','X','U'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Sage versus Python\n",
    "\n",
    "You should know that even if Sage is built on top of Python, these two\n",
    "languages do not behave exactly the same. To run a cell directly in\n",
    "Python start it with `%%python` as in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%python\n",
    "print(3^5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What is the answer to the above command in Sage?\n",
    "\n",
    "To make life of mathematician easier, Sage has a *preparser*: before a\n",
    "command is sent to Python it is transformed. You can call explicitely\n",
    "the preparser on a string to know what it does:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "preparse(\"3^5\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What does the preparser do to integers? To the operator `^`?\n",
    "\n",
    "You learnt that Python integers and Sage integers are not the same\n",
    "thing!\n",
    "\n",
    "## Some links\n",
    "\n",
    "### The essentials\n",
    "\n",
    "-   The main website: <http://www.sagemath.org/>\n",
    "-   A forum to ask your questions about Sage: <http://ask.sagemath.org>\n",
    "-   A book \"Calcul mathématique avec Sage\"/\"Computational Mathematics\n",
    "    with SageMath\"/\"Rechnen mit Sage\", a book about Sage (in french,\n",
    "    english and german): <http://sagebook.gforge.inria.fr/>\n",
    "\n",
    "### Introductory tutorials\n",
    "\n",
    "If you just start with Sage, it is a good idea to work on the 6\n",
    "Programming worksheets (\"First steps with Sage\", \"Learn about for\n",
    "loops\", etc) available on the [wiki](https://wiki.sagemath.org/days100).\n",
    "\n",
    "You can also have a look at the Sage documentation. These documents are\n",
    "part of Sage and you can access them from the Jupyter notebook by\n",
    "clicking on \"Help\" -\\> \"Thematic Tutorials\". They are also available at\n",
    "<http://doc.sagemath.org/>\n",
    "\n",
    "------------------------------------------------------------------------\n",
    "\n",
    "Authors  \n",
    "-   Thierry Monteil\n",
    "-   Vincent Delecroix\n",
    "\n",
    "License  \n",
    "CC BY-SA 3.0"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 8.9.beta3",
   "language": "sage",
   "name": "sagemath"
  },
  "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}