{
"nbformat_minor": 2,
"nbformat": 4,
"cells": [
{
"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",
"# Tutorial: Real and complex numbers in Sage\n",
"\n",
"Authors \n",
"- Thierry Monteil\n",
"\n",
"License \n",
"CC BY-SA 3.0\n",
"\n",
"Computers can not represent all real or complex numbers. With unlimited memory, you could represent all integers and rationals but not the real numbers since there are uncountably many of them. Sage proposes different solutions to make computations with real numbers. Each different way corresponds to a parent in the categories `Fields` or `Rings` though they do not necessarily satisfy their respective axiomatic. These different solutions have very different features and it is important to know to which kind of problem which representation to choose.\n",
"\n",
"You should know a bit about coercion before starting this tutorial, so it is advised to read the `parent_element_coercion.en.rst` worksheet first.\n",
"\n",
"## The big zoo\n",
"\n",
"Here is a rough qualitative summary of the main representations that will appear in this tutorial:\n",
"\n",
"| Kind | Name | Short | Complex equivalent | Short | Reliable | Fast | Expressive | Exact | Precision |\n",
"|-----------|---------------------------|-------|------------------------------|---------|----------|---------|------------|-------|-------------------------------|\n",
"| Numerical | `RealDoubleField()` | `RDF` | `ComplexDoubleField()` | `CDF` | `+++` | `+++++` | `+` | No | bounded by 53 bits |\n",
"| | `RealField(prec)` | `RR` | `ComplexField(prec)` | `CC` | `+++` | `+++` | `++` | No | finite unbounded non-adaptive |\n",
"| | `RealIntervalField(prec)` | `RIF` | `ComplexIntervalField(prec)` | `CIF` | `+++++` | `+++` | `++` | No | finite unbounded non-adaptive |\n",
"| | `RealBallField(prec)` | `RBF` | `ComplexBallField(prec)` | `CBF` | `+++++` | `+++` | `++` | No | finite unbounded non-adaptive |\n",
"| Algebraic | `AlgebraicRealField()` | `AA` | `AlgebraicField()` | `QQbar` | `+++++` | `+` | `++++` | Yes | infinite |\n",
"| | `NumberField(poly)` | none | itself | none | `+++++` | `++` | `++++` | Yes | infinite |\n",
"| | `QuadraticField(n)` | none | itself | none | `+++++` | `+++` | `+++` | Yes | infinite |\n",
"| | `RationalField()` | `QQ` | `QuadraticField(-1)` | none | `+++++` | `++++` | `+++` | Yes | infinite |\n",
"| Symbolic | `SymbolicRing()` | `SR` | itself | `SR` | `+` | `+` | `+++++` | No | depends on construction |\n",
"\n",
"Numerical representations are fast but inexact, they require some special care since the precision may become fuzzy along a computation.\n",
"\n",
"Algebraic representations are exact and reliable, but slower and able to represent only a class of numbers.\n",
"\n",
"Symbolic representations are very versatile, e.g. they can represent numbers such as $\\sqrt{log(\\pi)+sin(42)}$, but slow and sometimes not reliable.\n",
"\n",
"Let us inspect each kind of representation further.\n",
"\n",
"**TODO**\n",
"\n",
"Find a better word for \"kind\" that is not \"type\", \"group\", \"category\", \"class\" (perhaps \"taxon\", \"family\", \"sort\" ?).\n",
"\n",
"## Numerical representations\n",
"\n",
"### Floating-point numbers\n",
"\n",
"Floating-point numbers are real numbers of the form $(-1)^s 0.d_0 d_1 d_2 \\dots\n",
"d_{p-1} 2^{e}$ where:\n",
"\n",
"- $s \\in \\{0,1\\}$ is the *sign*,\n",
"- $d_i \\in \\{0,1\\}$ (with $d_0=1$) are the *digits*,\n",
"- the integer $d_0 d_1 d_2 \\dots d_{p-1}$ is the *mantissa*,\n",
"- $p$ is the *precision*,\n",
"- $e$ is the *exponent*.\n",
"\n",
"Floating-point numbers of given precision do not form a field in the mathematical sense since they are not stable under the classical operations. So when adding or multiplying two floating-point numbers, Sage has to round the result to return a floating-point number. In particular, the results are approximate and the error can grow along a computation.\n",
"\n",
"Some first warnings about floating-point arithmetics can be found on the guide: \n",
"\n",
"Error propagation is sometimes unavoidable, but it is possible to know a bound on it by using the following secure fields:\n",
"\n",
"There exists various floating-point implementations of numbers in Sage:\n",
"\n",
"#### Real Double Field\n",
"\n",
"Elements of `RDF = RealDoubleField()` are the floating-point numbers in double precision from the processor (those you use when you code in C, similar to Python $float$), see sage.rings.real\\_double. The computations using those numbers are very fast and most optimized libraries that work with floating-point numbers use this representation, so you will benefit from them if the number involved in your floating-point constructions (matrices, polynomials,...) belong to `RDF`. The mantissa has 53 bits of precision, the exponent is encoded on 11 bits from -1023 to 1024. Note that since the first bit of the mantissa is always $1$, we only have to store $1 + 52 + 11 = 64$ bits as expected. During the computations, the rounding is done toward the closest floating-point number.\n",
"\n",
"#### Real Fields\n",
"\n",
"Elements of `RealField(prec)` are the floating-point numbers with `prec` bits of precision, \"emulated\" by the MPFR library. By default, the mantissa has $53$ bits of precision (and `RealField(prec=53)` can be shortened as `RR`), the exponent is encoded on 32 bits (so that you will unlikely overflow), see sage.rings.real\\_mpfr. The computations are slower than in `RDF` and most algorithms on objects whose base ring is `RealField(prec)` will be pretty naive. However, it is possible to extend the precision as well as the type of rounding (towards $+\\infty$, towards $-\\infty$, towards $0$), for example `RealField(100, rnd='RNDZ')` is the \"field\" of real floating-point numbers with 100 bits of precision with rounding towards zero.\n",
"\n",
"> **note**\n",
">\n",
"> While sometimes needed, increasing the precision also increases the computation time.\n",
"\n",
"> **warning**\n",
">\n",
"> Note that despite its generic name, `RealField(prec)` or `RR` is *not* an abstraction of the field of real numbers, as can be seen by the following section.\n",
"\n",
"#### Unreal floats\n",
"\n",
"Besides not representing all real numbers, `RDF` and `RR` contain three additional unreal elements that allow the `(+,-,*,/)` operations to be always defined (`NaN` stands for *Not a Number*):"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RDF(0) / RDF(0)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"NaN\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RDF(1) / RDF(0)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"+infinity\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(-1) / RR(0)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"-infinity\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RDF(infinity) - RDF(infinity)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"NaN\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"NaN in RR and Infinity in RR and -Infinity in RR"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"#### Back and forth in precision\n",
"\n",
"Note that it is is general useless to improve the precision during a computation. This explains why the coercion is done toward the smallest precision:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"R = RealField(200)\n",
"S = RealField(100)\n",
"a = R(2).sqrt()\n",
"b = S(pi)\n",
"a+b == R(sqrt(2)+pi)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a+R(b) == R(sqrt(2)+pi)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"False"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"> **note**\n",
">\n",
"> For each element of those two fields, there exists a `.sign_mantissa_exponent()` method that returns the corresponding triple.\n",
"\n",
"### Floating-point numbers with a control on the imprecision\n",
"\n",
"#### Real Interval Fields\n",
"\n",
"Elements of `RealIntervalField(prec)` are pairs of elements of `RealField(prec)`, which should be interpreted as the two endpoints of an interval containing the represented real number. `RealIntervalField(53)` can be shortened as `RIF`.\n",
"\n",
"Operations defined over the `RealIntervalField` implementation must offer the *containment guaranty*: for every mathematical function $f : \\mathbb{R}^d \\to\n",
"\\mathbb{R}$ that admits an implementation `f_rif` over some `RealIntervalField(prec)`, we have the property that for every element `x = [(a_0, b_0), (a_1,b_1), ..., (a_(d-1),b_(d-1))]` in `RealIntervalField(prec)^d`, if `f_rif(x) = [a', b']`, then $f([a_0, b_0]\n",
"\\times [a_1,b_1] \\times \\dots \\times [a_(d-1),b_(d-1)])\\subseteq [a', b']$. In particular, the left endpoint of the interval is rounded towards $-\\infty$ and the right endpoint is rounded towards $+\\infty$\n",
"\n",
"#### Real Ball Fields\n",
"\n",
"Elements of `RealBallField(prec)` are pairs `(c,r)`, where `c` is a floating-point number with `prec` bits of precision, and where `r` is a floating-point number with 30 bits of precision. The high precision `c` and the low precision number `r` should be understood as the center and the radius (respectively) of a ball containing the represented real number. `RealBallField(53)` can be shortened as `RBF`.\n",
"\n",
"Operations defined over the `RealBallField` implementation must offer a similar containment guaranty as for `RealIntervalField` (replace intervals with balls).\n",
"\n",
"#### Comparison between both representations\n",
"\n",
"The main advantage of the ball implementation over interval arithmetics is that, when the precision is high, the error is still encoded with a low precision (since only the order of magnitude matters for the error bounds), so that the loss of speed compared to non-guaranteed floating-point arithmetics with the same precision is about $1+\\varepsilon$ (and not $2$ as in the case of interval arithmetics, where two precise numbers have to be computed).\n",
"\n",
"Also, it is closer to the mathematical $\\varepsilon-\\delta$ definitions and distance-based error estimations.\n",
"\n",
"However, for large intervals and for evaluating monotonic functions, interval arithmetics is usually more precise and easier to implement (just round the images of the endpoints towards the right direction). It is also well-adapted for implementing dichotomy-based precision refinements.\n",
"\n",
"### Which numerical representation to choose\n",
"\n",
"Sage interfaces various libraries written in C, C++, Fortran, etc. which use the processor double precision floats. Hence, the use of `RDF` has the advantage of being faster and allowing the use of those libraries.\n",
"\n",
"If you use other floating-point representations, Sage will mostly use generic algorithms, which are usually slower and numerically less stable. However, those representations allow you to express real numbers with a higher precision.\n",
"\n",
"Regarding speed, `RealBallField(prec)` and `RealField(prec)` are similar, while `RealIntervalField(prec)` is twice slower. However, the `RealBallField(prec)` is a new feature of Sage, hence some methods might raise a `NotImplementedError`.\n",
"\n",
"**Summary**\n",
"\n",
"So, if you need speed and the power of Sage libraries: use `RDF`.\n",
"\n",
"If you need more precision, start with `RealBallField(prec)`, and try `RealIntervalField(prec)` in case the imprecision becomes too large (alternatively: increase the precision or try to think about the order and the expanding nature of your operations). Use `RealField(prec)` if the methods you wanted to use with `RealBallField(prec)` are not implemented.\n",
"\n",
"### Exercises\n",
"\n",
"#### Binomial\n",
"\n",
"From the three following values:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 1.0\n",
"b = 10.0^4\n",
"c = 1.0"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"solve the $ax^2+bx+c = 0$ equation by the method you learned in highschool involving discriminant, and compare the sum and product of the solution you found to their theoretical values $-b/a$ and $c/a$."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Same question by letting Sage solve this equation for you."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Even on such a simple computation, you have been victim of a phenomenon of *error amplification*. Another striking example is the following:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 10000.0\n",
"b = -9999.5\n",
"c = 0.1\n",
"a+c+b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.600000000000364"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a+b+c"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.600000000000000"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Note that, by construction and following the IEEE 754 specification, the `+` (and `*`) law is commutative in `RDF` and `RR` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"all([a+b==b+a, a+c==c+a, b+c==c+b, (a+b)+c==c+(a+b), (a+c)+b==b+(a+c), (b+c)+a==a+(b+c)])"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Actually the associativity is not satisfied here:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a+b)+c == a+(b+c)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"False"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Try to find, by random sampling, 3 elements $a,b,c$ in `RDF` such that `(a+b)+c != a+(b+c)` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Why is it so hard ? What makes the previous example work ?\n",
"\n",
"#### Interval arithmetics\n",
"\n",
"\"Prove\" using `RealIntervalField` that $cos(2\\pi)=1$ with a high precision (say $10^{-100}$)."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here "
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Explain the following behaviour:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"sqrt(0.1)^2 == 0.1"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"sqrt(RIF(0.1))^2 == RIF(0.1)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"False"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Provide a relation between `sqrt(RIF(0.1))^2` and `RIF(0.1)`"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Note that the ordering of the computations has an influence on the result, but also on its precision. With the previous example, we have:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = RIF(10000.0)\n",
"b = RIF(9999.5)\n",
"c = RIF(0.1)\n",
"a+c-b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.60000000000?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a-b+c"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.6000000000000000?\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a+c-b).endpoints()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(0.599999999998544, 0.600000000000364)"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a+c-b).absolute_diameter()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.81898940354586e-12\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a-b+c).endpoints()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(0.599999999999999, 0.600000000000001)"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a-b+c).absolute_diameter()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.11022302462516e-16"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Before evaluating them, could you predict and explain the behaviour of the following commands ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# 1/2 in RIF"
],
"outputs": [],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# 1/3 in RIF"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"*Hint*: you can look at the documentation of the `.__contains__()` method of `RIF`.\n",
"\n",
"A Sage developer wants to write the real function $f$ that maps $0.01$ to $1$ and that is the identity elsewhere, and suggests the following implementation:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"def f_unstable(x):\n",
" if x == 0.01:\n",
" return 1\n",
" else:\n",
" return x"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Prove that this implementation does not offers the containment guaranty when used on the `RealIntervalField`.\n",
"\n",
"Could you help this lost developer to find a correct implementation so that it could be included in Sage source code ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"def f_correct(x):\n",
" # edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### Divergent series ?\n",
"\n",
"The series $u_n := 1 + 1/2 + 1/3 + 1/4 + \\dots + 1/n = \\sum_{k=1}^{n} 1/k$ is known to diverge on $\\mathbb{R}$.\n",
"\n",
"Explore its convergence on some $Real Fields$ (start with small precision). Explore the difference with the (mathematically equal) series $v_n := 1/n +\n",
"1/(n-1) + ... + 1/2 + 1$.\n",
"\n",
"Explain, formulate conjectures and try to prove them."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### A picture\n",
"\n",
"Draw (by hand) a picture of the elements of `RealField(3)`.\n",
"\n",
"#### Convergence of a Markov Chain\n",
"\n",
"Let $M$ be the rational bistochastic matrix defined by:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"M = matrix(QQ, [[1/6, 1/2, 1/3], [1/2, 1/4, 1/4], [1/3, 1/4, 5/12]])\n",
"M"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"[ 1/6 1/2 1/3]\n",
"[ 1/2 1/4 1/4]\n",
"[ 1/3 1/4 5/12]"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Could you guess the limit of $M^n$ when $n$ goes to $+\\infty$ ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Could you prove that this limit exists and give its exact value ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### Repulsive and attractive fixed points\n",
"\n",
"Define the maps $f:x\\to 4x-1$ and $g:x\\to x/3+1$."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Let $u_n$ and $v_n$ be the sequences defined recursively by $u_0 = RDF(1/3)$, $u_{n+1} = f(u_n)$ and $v_0 = RDF(3/2)$, $v_{n+1} = g(v_n)$.\n",
"\n",
"Compute the first 100 values of the sequences $u_n$ and $v_n$."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"For which reason you should not worry to find the fixed point of the function $h:x\\to sin(x + 1/2)$ ? Compute the value of this fixed point with 50 bits of precision."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### Decimals of pi\n",
"\n",
"What is the 150th decimal of $\\pi$ ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### Money, money, money"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"one_cent = 0.01\n",
"my_money = 0\n",
"for i in range(100):\n",
" my_money += one_cent\n",
"my_money > 1"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"What is the origin of this behaviour ?\n",
"\n",
"How to fix the issue in Sage (respecively un pure Python, compare both) ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Which risk does computer science take if finance and high-frequency trading stays on power for a while ?\n",
"\n",
"Provide the first 30 bits of the binary expansion of $0.01$."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"#### Huge matrix computation\n",
"\n",
"Construct a dense matrix of size `500*500` with real coefficients randomly distributed in $[0,1]$."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Take the 100th power of this matrix in less than one second."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Comment.\n",
"\n",
"#### Exponents"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RDF(1/2^10000) == 0"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(1/2^10000) == 0"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"False"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Comment.\n",
"\n",
"#### A strange impression"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(1/pi)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.318309886183791"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RDF(1/pi)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.3183098861837907"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Before evaluating the next command, could you tell what will be the answer of the following command ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# RR(1/pi) == RDF(1/pi)"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"What is the cause of the behaviour of the first two commands ?\n",
"\n",
"#### Rounding\n",
"\n",
"Let us define:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"n = 5.0000000000000001"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Explain each of the following behaviours:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"n"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"5.000000000000000"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"n.str(truncate=False)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"'5.000000000000000111'"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"5.000000000000000111 == 5.0000000000000001"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"## Algebraic (exact) representations\n",
"\n",
"Some subfields of the real numbers can be represented in an exact way in Sage. They are pretty safe (no rounding problem, well defined semantics) through sometimes slow.\n",
"\n",
"As a general rule, the smallest the field is (for the inclusion), the fastest it is, which means that, in principle, working within rational field is faster than within quadratic fields, which is faster than within cyclotomic fields, which is faster than within number fields, which is faster that within the universal cyclotomic field, that is faster than within the whole algebraic field.\n",
"\n",
"Here is a short overview:\n",
"\n",
"### Rational numbers\n",
"\n",
"`QQ = RationalField()` stands for the field of rational numbers. Computations are exact and only limited by the available memory. You should note that the numerators/denominators could become very large along a computation, which might slow down this latter."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(123/321)^10000"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"It is possible to \"approximate\" a floating-point real number by rational numbers with the methods `nearby_rational()`, `simplest_rational()`, `exact_rational()`."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"If `x` belongs to `RR`, what is the result of `QQ(x)` ? Is it a coercion ?\n",
"\n",
"See the page \n",
"\n",
"### Algebraic number field\n",
"\n",
"`AA = AlgebraicRealField()` stands for te field of real algebraic numbers."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = AA(2).sqrt() + AA(7).sqrt()\n",
"a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.059964873437685?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Here, the number `a` is *printed* as its approximated numerical value, but the question mark at the end indicates that it is not internally a floating-point number.\n",
"\n",
"Basically, an algebraic number is represented as the root of an integer polynomial, and a small interval surrounding the appropriate root:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.minpoly()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"x^4 - 18*x^2 + 25\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.n(digits=50)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.0599648734376856393033044778489585042799310584594"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.interval_diameter(1e-40)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.059964873437685639303304477848958504279931058459398253545014197191801301693?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Computing such a representation is costly, hence Sage is lazy and computes it only when needed:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b = QQbar(2).sqrt()\n",
"c = b^2\n",
"c"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"2.000000000000000?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Here, Sage did not compute the minimal polynomial nor an enclosing interval of `c`, it just knows `c` as \"the square of `b`\", in particular, it does not know yet that `c` is the integer `2`. However, such an \"exactification\" is computed when required (and kept internally):"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"c.is_integer()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"c"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"2"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"You can force the exactification with the `exactify` method.\n",
"\n",
"For more details, see the page \n",
"\n",
"### Number fields\n",
"\n",
"See the pages:\n",
"\n",
"- \n",
"- \n",
"- \n",
"\n",
"> **warning**\n",
">\n",
"> You should be aware that number field are purely algebraic objects that are *not* embedded into the complex plane by default. Do to this, you have to define the embedding you want to work with (if you need one).\n",
"\n",
"Note that the quadratic number fields have a dedicated faster implementation.\n",
"\n",
"What is the difference between those two objects ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"q = QuadraticField(-1,'I')\n",
"qq = QQ[I]"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"### Cyclotomic fields (and the universal cyclotomic field)\n",
"\n",
"See the pages:\n",
"\n",
"- \n",
"- \n",
"\n",
"## Symbolic representations\n",
"\n",
"`SR = sage.symbolic.ring.SymbolicRing()` is the set of symbolic expressions.\n",
"\n",
"`SR` is pretty fuzzy but quite handy (here expressivity is prefered over standardization and the number of real numbers representable within `SR` is constantly increasing). This \"ring\" contains most mathematical constants such as $e$, $\\pi$, $cos(12)$, $\\sqrt{2}$, $\\gamma$ (Euler constant), ... and allows some computations such as:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"e^(I*pi)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"-1"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"You should not consider `SR` as reliable:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"bool(567663097408 * pi - 1783366216531 == 0)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"> **warning**\n",
">\n",
"> `SR` also contains symbolic expressions that are not real or complex numbers:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"e = cos(x)^2\n",
"e.derivative()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"-2*cos(x)*sin(x)"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"_ in SR"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"`SR` represents numbers as symbolic expressions, whose length (as a formula) can grow along a computation. This phenomenon can slow down the computations. If our aim is to end with a numerical value, you should use a floating-point representation, as they keep a constant length."
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = exp(2)\n",
"f = lambda x : sqrt(x+1)^x\n",
"for i in range(10): \n",
" a = f(a)\n",
"a"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Of course, using floating-point arithmetics will not necessarily solve all your problems since you will have to ensure the numerical stability of such iterated computation (see the sections above).\n",
"\n",
"`SR` always answer, even when it does not know the answer, for example:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# missing example"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"Hence, computations involving the symbolic ring should not be considered as reliable. In particular, `SR` is a dead-end for the coercion among numbers (and can add pears and apples formally):"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"0.1 + pi + QQbar(2).sqrt()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"pi + 1.51421356237310"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"_.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Symbolic Ring"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Here are some issues that can appear with the symbolic ring, when you can represent your numbers in a different way, we strongly encourage you to do so:\n",
"\n",
"- \n",
"- \n",
"- \n",
"- \n",
"- \n",
"- \n",
"\n",
"## Conversions, bridges, exactification\n",
"\n",
"**TODO**\n",
"\n",
"Should this be a separate section, or should it appear in the rational/algebraic sections ?\n",
"\n",
"Along a computation, coercion goes towards the less precise numbers:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 1/10 ; a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1/10"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Rational Field"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b = QQbar(2).sqrt() ; b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.414213562373095?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Algebraic Field"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"c = a+b ; c"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.514213562373096?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"c.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Algebraic Field"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"d = RDF(0.1) ; d"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.1"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"d.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Real Double Field"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"e = c+d ; e"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.614213562373095"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"e.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Complex Double Field"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"f = pi ; f"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"pi"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"f.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Symbolic Ring"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"g = e+f ; g"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"pi + 1.614213562373095"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"g.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Symbolic Ring"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"But you might want to swim upstream. Given a floating point approximation of a real number, there are various ways to recover an exact representation it *may* come from.\n",
"\n",
"### From numerical to rational\n",
"\n",
"An element `x` of `RealField(prec)` can be transformed into a rational in different ways:\n",
"\n",
"- `x.exact_rational()` returns the exact rational representation of the floating-point number `x`. In particular its denominator is a power of\n",
" 1. It contains the same information as `.sign_mantissa_exponent` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = RR(pi) ; a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"3.14159265358979"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.sign_mantissa_exponent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(1, 7074237752028440, -51)"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"s,m,e = _\n",
"s*m*2^(e)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"884279719003555/281474976710656"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b = a.exact_rational() ; b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"884279719003555/281474976710656"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"- `x.simplest_rational()` returns the simplest rational (i.e. with smallest denominator) that will become equal to `x` when converted into the parent of `x` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = RR(pi) ; a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"3.14159265358979"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b = a.simplest_rational() ; b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"245850922/78256779"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a == b"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(b)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"3.14159265358979"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"- `nearby_rational(max_error=e)` finds the simplest rational that is at distance at most `e` of `x` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(pi).nearby_rational(max_error=0.01)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"22/7"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"- `nearby_rational(max_denominator=d)` finds the closest rational to `x` with denominator at most `d` :"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RR(pi).nearby_rational(max_denominator=120)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"355/113"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"### From numerical to algebraic\n",
"\n",
"If `x` is an element of `RealField(prec)`, `x.algebraic_dependency(d)` returns an integer polynomial `P` of degree at most `d` such that `P(x)` is almost zero:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 3.0.nth_root(3)\n",
"a += 5*a.ulp()\n",
"a.algebraic_dependency(3)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"x^3 - 3"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"### From algebraic field to (smaller, faster) embedded number fields\n",
"\n",
"Algebraic numbers have an `as_number_field_element` method:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = QQbar(2).sqrt() + QQbar(3).nth_root(3)\n",
"a.as_number_field_element()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(Number Field in a with defining polynomial y^6 - 6*y^4 - 6*y^3 + 12*y^2 - 36*y + 1,\n",
" 96/755*a^5 - 54/755*a^4 + 128/151*a^3 + 936/755*a^2 - 1003/755*a + 2184/755,\n",
" Ring morphism:\n",
" From: Number Field in a with defining polynomial y^6 - 6*y^4 - 6*y^3 + 12*y^2 - 36*y + 1\n",
" To: Algebraic Real Field\n",
" Defn: a |--> -0.02803600793431334?)\n"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"field, number, morphism = _\n",
"morphism(number) == a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"For more than one algebraic number, the previous method generalizes with the `number_field_elements_from_algebraics` function:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"from sage.rings.qqbar import number_field_elements_from_algebraics\n",
"sqrt2 = AA(2).sqrt()\n",
"sqrt3 = AA(3).sqrt()\n",
"number_field_elements_from_algebraics([sqrt2, sqrt3])"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(Number Field in a with defining polynomial y^4 - 4*y^2 + 1,\n",
" [-a^3 + 3*a, -a^2 + 2],\n",
" Ring morphism:\n",
" From: Number Field in a with defining polynomial y^4 - 4*y^2 + 1\n",
" To: Algebraic Real Field\n",
" Defn: a |--> 0.5176380902050415?)"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"field, numbers, morphism = _\n",
"morphism(numbers[0]) == sqrt2"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"### From algebraic to symbolic\n",
"\n",
"Algebraic numbers have a `radical_expression` method:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"x = polygen(QQ,'x')\n",
"P = x^3-x^2-x-1\n",
"P"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"x^3 - x^2 - x - 1"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"P.roots(QQbar)[0][0]"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.839286755214161?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"_.radical_expression()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"(1/9*sqrt(11)*sqrt(3) + 19/27)^(1/3) + 4/9/(1/9*sqrt(11)*sqrt(3) + 19/27)^(1/3) + 1/3"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"> **note**\n",
">\n",
"> Note that the `radical_expression` method does not implement the whole Galois theory yet, meaning that some algebraic numbers that can be represented by radicals, might not be modified that way:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = QQbar(2).sqrt() + QQbar(3).nth_root(3)\n",
"a.radical_expression()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"2.856463132680504?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"### Exercise: reverse engineering\n",
"\n",
"Suppose that, with a numerical experiment, you have approached an irrational algebraic number by:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 3.14626436994197234"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"What is this algebraic number likely to be ?"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"# edit here"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"## Complex numbers\n",
"\n",
"Most real number representation have their complex equivalent.\n",
"\n",
"- `RR` becomes `CC=ComplexField()`\n",
"- `RDF` becomes `CDF=ComplexDoubleField()`\n",
"- `RIF` becomes `CIF=ComplexIntervalField()`\n",
"- `RBF` becomes `CBF=ComplexBallField()`\n",
"- `QQ` becomes `QQ[I]` (Gauss integers)\n",
"- `AA` becomes `QQbar` (the algebraic closure of `QQ`)\n",
"- `SR` remains `SR` (which contains much more than just complex numbers)\n",
"\n",
"Note that by default, `I` and `i` belong to the `Symbolic Ring`, so you will have to convert it, e.g. `QQbar(I)`, `CDF(I)`, ... or obtain it as `QQbar.gen()`, `QQbar.gen()`, ... . If you overwrite it, you can recover it with:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"from sage.symbolic.all import i as I"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"### Diamonds are not round\n",
"\n",
"We saw previously that contracting makes computations safe. Consider the following example:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = CBF(1/10+I*4/5)\n",
"b = CBF((1+I)*2/3)\n",
"a.diameter()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.4408921e-16"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b.diameter()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.4408921e-16"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"b.abs() < 1"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"True"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"(a*b).diameter()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"1.1768364e-15"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Explain this loss of precision.\n",
"\n",
"*Hint*: draw a picture !\n",
"\n",
"Complex ball elements are pairs of real ball elements, see \n",
"\n",
"## Transitional representations (under the hood)\n",
"\n",
"Some \"transitional representations\" were not included in the big zoo for simplicity, and because they exist mostly as wrappers for internal stuff.\n",
"\n",
"### Real literals: between your fingers and Sage\n",
"\n",
"When the user writes:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = 0.1"
],
"outputs": [],
"metadata": {}
},
{
"source": [
"this defines a floating-point number though the precision is not explicitely set.\n",
"\n",
"Sage defines some default precision (depending on the number of digits that are provided, 53 in the current case):"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"type(a)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
""
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a.parent()"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"Real Field with 53 bits of precision"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"but it allows casting into higher precision rings, without losing the precision from the default choice:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RealField(200)(0.1)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.10000000000000000000000000000000000000000000000000000000000"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RealField(200)(RDF(0.1))"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"0.10000000000000000555111512312578270211815834045410156250000"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"### Real Lazy Field: between exact and numerical representations\n",
"\n",
"Elements of `RLF = RealLazyField()` behaves like floating-point numbers, This represen\n",
"\n",
"while\n",
"\n",
"`RLF = RealLazyField()` is lazy, in the sense that it doesn't really do anything but simply sits between exact rings of characteristic 0 and the real numbers. The values are actually computed when they are cast into a field of fixed precision.\n",
"\n",
"`RLF = RealLazyField()` like `RealField()` but can wait that the user fixes the precision to get one. "
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a = RLF(pi) + 2\n",
"a"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"5.141592653589794?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RealField(100)(a)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"5.1415926535897932384626433833"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RealField(150)(a)"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"5.1415926535897932384626433832795028841971694"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"The precision of `a` can be increased on demand since, internally, `a` keeps the exact representation it comes from:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"a._value"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"pi + 2"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"source": [
"Its use is sometimes tricky, and shares some weirdnesses of the symbolic ring since it also warps its elements:"
],
"cell_type": "markdown",
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"RLF(pi) + RLF(QQbar(2).sqrt())"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"4.555806215962889?"
]
},
"metadata": {}
}
],
"metadata": {}
},
{
"execution_count": null,
"cell_type": "code",
"source": [
"_._value"
],
"outputs": [
{
"execution_count": 1,
"output_type": "execute_result",
"data": {
"text/plain": [
"pi + 1.414213562373095?"
]
},
"metadata": {}
}
],
"metadata": {}
}
],
"metadata": {
"kernelspec": {
"display_name": "sagemath",
"name": "sagemath"
}
}
}