{ "cells": [ { "cell_type": "markdown", "id": "70df0d8a-df83-4330-8bcf-0d8fd63aac59", "metadata": {}, "source": [ "Submit your answers to gradescope by 3:40pm." ] }, { "cell_type": "markdown", "id": "26038e7a-565c-4ddb-94a2-78b271eb1c0e", "metadata": {}, "source": [ "# question 1" ] }, { "cell_type": "markdown", "id": "e28f8844-581f-4567-94d1-4f86aae8da1b", "metadata": {}, "source": [ "Compute the first 10 rational approximations of\n", "$\\displaystyle \\cos(\\pi/5)$, with increasing accuracy,\n", "starting at 1 decimal place of accuracy,\n", "going up to 10 decimal places.\n", "\n", "Verify the accuracy of the approximations via their relative errors." ] }, { "cell_type": "markdown", "id": "3eaa057f-ad15-4504-8178-e1e425e5e759", "metadata": {}, "source": [ "## answer to question 1" ] }, { "cell_type": "code", "execution_count": 1, "id": "9a000844-15d2-4f7d-a9d7-f4794d720548", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\frac{1}{4} \\, \\sqrt{5} + \\frac{1}{4}\\)" ], "text/latex": [ "$\\displaystyle \\frac{1}{4} \\, \\sqrt{5} + \\frac{1}{4}$" ], "text/plain": [ "1/4*sqrt(5) + 1/4" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "cpi5 = cos(pi/5)\n", "show(cpi5)" ] }, { "cell_type": "code", "execution_count": 2, "id": "930413bc-140c-4832-93d0-d7b1f5641a8d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[13/16,\n", " 38/47,\n", " 161/199,\n", " 377/466,\n", " 1292/1597,\n", " 4181/5168,\n", " 17711/21892,\n", " 23184/28657,\n", " 98209/121393,\n", " 416020/514229]" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A = [QQ(cpi5.n(digits=k)) for k in range(1, 11)]\n", "A" ] }, { "cell_type": "markdown", "id": "eb9f9293-f868-41ae-aff4-b69e9352d496", "metadata": {}, "source": [ "To verify the accuracy, we compute the relative errors." ] }, { "cell_type": "code", "execution_count": 3, "id": "d4f1f298-3a2f-4720-9b8b-734e8cebaf2e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[0.00430523171857909,\n", " 0.000625890532085027,\n", " 0.0000348963691764936,\n", " 5.09116433322903e-6,\n", " 1.08372117528138e-7,\n", " 4.13944698704118e-8,\n", " 2.30683469700907e-9,\n", " 3.36562735072396e-10,\n", " 1.87559251193299e-11,\n", " 1.04528938667147e-12]" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "errs = [abs(RR(A[k] - cpi5))/RR(cpi5) for k in range(10)]\n", "errs" ] }, { "cell_type": "markdown", "id": "09441c6b-d458-41f0-89f9-0928e797dae0", "metadata": {}, "source": [ "For a better representation of the errors, we convert to the scientific notation." ] }, { "cell_type": "code", "execution_count": 4, "id": "db1a807c-d826-4db4-9347-8f3d01813f15", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "['4.31E-3',\n", " '6.26E-4',\n", " '3.49E-5',\n", " '5.09E-6',\n", " '1.08E-7',\n", " '4.14E-8',\n", " '2.31E-9',\n", " '3.37E-10',\n", " '1.88E-11',\n", " '1.05E-12']" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "[f\"{error:.2E}\" for error in errs]" ] }, { "cell_type": "markdown", "id": "4b9a34b8-2766-4edc-9341-4722756f6282", "metadata": {}, "source": [ "The relative errors match the required accuracy of the rational approximations." ] }, { "cell_type": "markdown", "id": "d665a5ec-6ac6-4b80-a72d-fa3f7b1ffbde", "metadata": {}, "source": [ "# question 2" ] }, { "cell_type": "markdown", "id": "b25f93f0-1f56-44e0-b9f8-162303ff6630", "metadata": {}, "source": [ "Compute the square roots of the complex number $c = 7 + 3 I$.\n", "\n", "For each of the square roots $r$ verify symbolically that $r^2 = c$." ] }, { "cell_type": "markdown", "id": "16b84ef8-7764-4bf0-bf7e-c5250dd678ec", "metadata": {}, "source": [ "## answer to question 2" ] }, { "cell_type": "code", "execution_count": 5, "id": "e7da5aeb-f400-43ec-8b2a-43c2d486bb60", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 6, "id": "ef8517fd-9765-4756-8726-6d09ab65d9fe", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle 3 i + 7\\)" ], "text/latex": [ "$\\displaystyle 3 i + 7$" ], "text/plain": [ "3*I + 7" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "c = 7 + 3*I\n", "show(c)" ] }, { "cell_type": "code", "execution_count": 7, "id": "32c14b6d-4612-49b0-87b6-2e49bb577a5c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[sqrt(3*I + 7), -sqrt(3*I + 7)]" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "r = sqrt(c, all=True)\n", "r" ] }, { "cell_type": "code", "execution_count": 8, "id": "922a8a84-24fa-4c16-80ff-b81ecd55048a", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[3*I + 7, 3*I + 7]" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "[a^2 for a in r]" ] }, { "cell_type": "markdown", "id": "2dc00a02-0bee-4bea-a930-a30f2de69350", "metadata": {}, "source": [ "Indeed, we obtain $c$ for all squares of the roots." ] }, { "cell_type": "markdown", "id": "0ea7dfe4-b2c1-4c8c-aa33-8946c2a45b3a", "metadata": {}, "source": [ "# question 3" ] }, { "cell_type": "markdown", "id": "1d1be447-a858-4edd-bf9d-5c308bb5f26a", "metadata": {}, "source": [ "Let $p = x^3 + 2 x + 5$ be a polynomial with coefficients\n", "in the integer number field modulo 11. \n", "Is $p$ irreducible over ${\\mathbb Z}_{11}$?\n", "\n", "Extend the coefficient field with a formal root of $p$\n", "so that, over this extended field, \n", "$p$ factors into a product of linear polynomials." ] }, { "cell_type": "markdown", "id": "91ca6d19-0182-4bf0-a63d-fbf0fd5d1843", "metadata": {}, "source": [ "## answer to question 3" ] }, { "cell_type": "code", "execution_count": 9, "id": "13ccb1b6-ac00-4558-9ce4-9762fd062d72", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 10, "id": "0cd37896-92ed-4a99-a460-b02016e0c0b2", "metadata": {}, "outputs": [], "source": [ "x = polygen(GF(11))" ] }, { "cell_type": "code", "execution_count": 11, "id": "3d3bbbc2-bd70-4f31-b752-fc264f9ea130", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle x^{3} + 2 x + 5\\)" ], "text/latex": [ "$\\displaystyle x^{3} + 2 x + 5$" ], "text/plain": [ "x^3 + 2*x + 5" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "p = x^3 + 2*x + 5\n", "show(p)" ] }, { "cell_type": "code", "execution_count": 12, "id": "46f624b3-2e73-476c-b4c5-1ae9a5dea5bd", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "p.is_irreducible()" ] }, { "cell_type": "markdown", "id": "140d4b8d-587d-4921-b9a1-9e190fcdbff4", "metadata": {}, "source": [ "As $p$ is not irreducible, we cannot use $p$ directly to extend the field of 11 elements. We will select the quadratic factor of $p$ to extend the field with a formal root." ] }, { "cell_type": "code", "execution_count": 13, "id": "5eab512f-d1d4-4438-b2fd-19ba5b19533e", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle (x + 7) \\cdot (x^{2} + 4 x + 7)\\)" ], "text/latex": [ "$\\displaystyle (x + 7) \\cdot (x^{2} + 4 x + 7)$" ], "text/plain": [ "(x + 7) * (x^2 + 4*x + 7)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "f = factor(p)\n", "show(f)" ] }, { "cell_type": "code", "execution_count": 14, "id": "28d5fe23-8228-437a-91ed-64d4c26795e0", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[(x + 7, 1), (x^2 + 4*x + 7, 1)]" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "L = list(f)\n", "L" ] }, { "cell_type": "code", "execution_count": 15, "id": "5211f931-4067-47ee-bb18-c87aea61fcd4", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle x^{2} + 4 x + 7\\)" ], "text/latex": [ "$\\displaystyle x^{2} + 4 x + 7$" ], "text/plain": [ "x^2 + 4*x + 7" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "f2 = L[1][0]\n", "show(f2)" ] }, { "cell_type": "code", "execution_count": 16, "id": "d8f91e44-9883-4598-b980-e00d82e9d38b", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Finite Field in a of size 11^2" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E. = GF(11).extension(f2)\n", "E" ] }, { "cell_type": "code", "execution_count": 17, "id": "28a5f6b2-5107-443d-9225-c303d02a3b27", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(x + 7) * (x + a + 4) * (x + 10*a)" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "factor(E[x](p))" ] }, { "cell_type": "markdown", "id": "04a1008b-76f9-4f36-a95d-34240c0e0004", "metadata": {}, "source": [ "Casting $p$ into the extended number field makes that it can be written as a product of three linear factors." ] }, { "cell_type": "markdown", "id": "90e8ed25-9991-4f4e-900e-c8ef57bce6c4", "metadata": {}, "source": [ "# question 4" ] }, { "cell_type": "markdown", "id": "8dcaecaa-d233-4c29-8e83-1046991d1b2a", "metadata": {}, "source": [ "What is delayed evaluation? \n", "\n", "Give a good application of a delayed evaluation in SageMath." ] }, { "cell_type": "markdown", "id": "79fc1f7d-c30c-4f33-bc4c-8fbc02908271", "metadata": {}, "source": [ "## answer to question 4" ] }, { "cell_type": "markdown", "id": "71e89b36-70c4-4cc6-a7f3-e28e92bc7ea6", "metadata": {}, "source": [ "Delayed evaluation postpones the evaluation of an expression, by placing a hold on the function call." ] }, { "cell_type": "code", "execution_count": 18, "id": "cc2a2d94-4bcb-43f5-846e-7ec9d75afa9e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1/4*sqrt(5) + 1/4" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cos(pi/5)" ] }, { "cell_type": "markdown", "id": "0cb2ca45-3ce6-49a8-85c8-3a17165573eb", "metadata": {}, "source": [ "We see that ``cos(pi/5)`` automatically evaluates to an expression involving the square root of 5. To prevent this from happening, we can place a hold on the function call, as is done below." ] }, { "cell_type": "code", "execution_count": 19, "id": "3cff5609-5957-42de-9d34-4182846fabd3", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\cos\\left(\\frac{1}{5} \\, \\pi\\right)\\)" ], "text/latex": [ "$\\displaystyle \\cos\\left(\\frac{1}{5} \\, \\pi\\right)$" ], "text/plain": [ "cos(1/5*pi)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "e = cos(pi/5, hold=True)\n", "show(e)" ] }, { "cell_type": "markdown", "id": "e6ccd743-4097-43f9-aea9-35af5b544bd9", "metadata": {}, "source": [ "Removing the hold evaluates the expression." ] }, { "cell_type": "code", "execution_count": 20, "id": "02db2020-052c-482b-9a4a-3cfe1277ea28", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1/4*sqrt(5) + 1/4" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "e.unhold()" ] }, { "cell_type": "markdown", "id": "b0ff39c0-b84f-4191-a7f5-434e0a309795", "metadata": {}, "source": [ "# question 5" ] }, { "cell_type": "markdown", "id": "1ab45241-0d93-4a69-8b5e-a13e50b52bc2", "metadata": {}, "source": [ "Draw the expression tree for\n", "$\\displaystyle \\frac{\\sin(x) + \\cos(x y) + 2}{x + y - 1}$." ] }, { "cell_type": "markdown", "id": "12fc58fa-1f7a-481a-b329-8b48776ed16f", "metadata": {}, "source": [ "## answer to question 5" ] }, { "cell_type": "code", "execution_count": 21, "id": "a35d90fe-bb64-4a8e-ba7f-a6942ed8077c", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 22, "id": "ee13a42c-b066-4ff3-9eaf-bbf37142d10e", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\frac{\\cos\\left(x y\\right) + \\sin\\left(x\\right) + 2}{x + y - 1}\\)" ], "text/latex": [ "$\\displaystyle \\frac{\\cos\\left(x y\\right) + \\sin\\left(x\\right) + 2}{x + y - 1}$" ], "text/plain": [ "(cos(x*y) + sin(x) + 2)/(x + y - 1)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x, y = var('x, y')\n", "expression = (sin(x) + cos(x*y) + 2)/(x + y - 1)\n", "show(expression)" ] }, { "cell_type": "code", "execution_count": 23, "id": "e18448ff-84d1-49bc-b9fb-6ff650ffa79f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[1/(x + y - 1), cos(x*y) + sin(x) + 2]" ] }, "execution_count": 23, "metadata": {}, "output_type": "execute_result" } ], "source": [ "e1 = expression.operands()\n", "e1" ] }, { "cell_type": "code", "execution_count": 24, "id": "af5653bf-879d-464c-96bb-7db0f41d2c2c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "expression.operator()" ] }, { "cell_type": "code", "execution_count": 25, "id": "81b0b404-363f-4630-b649-d03edb19ae4f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[x + y - 1, -1]" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "e1[0].operands()" ] }, { "cell_type": "code", "execution_count": 26, "id": "3df3308d-b099-40aa-b959-8c13dd25ba26", "metadata": {}, "outputs": [ { "data": { "text/plain": [ " _+__\n", " / / / \n", "x y -1" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Lx = LabelledOrderedTree([], label='x')\n", "Ly = LabelledOrderedTree([], label='y')\n", "Lm1 = LabelledOrderedTree([], label='-1')\n", "Ly = LabelledOrderedTree([], label='y')\n", "den = LabelledOrderedTree([Lx, Ly, Lm1], label='+')\n", "ascii_art(den)" ] }, { "cell_type": "code", "execution_count": 27, "id": "387eb1b6-6112-44b4-9af8-dd77f3c999c9", "metadata": {}, "outputs": [ { "data": { "text/plain": [ " _^___\n", " / / \n", " _+__ -1\n", " / / / \n", "x y -1 " ] }, "execution_count": 27, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f1 = LabelledOrderedTree([den, Lm1], label='^')\n", "ascii_art(f1)" ] }, { "cell_type": "code", "execution_count": 28, "id": "e585c57c-a640-411a-951b-192b0b7eb479", "metadata": {}, "outputs": [ { "data": { "text/plain": [ " cos\n", " |\n", " *_\n", " / /\n", "x y" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "xmy = LabelledOrderedTree([Lx, Ly], label='*')\n", "cosxy = LabelledOrderedTree([xmy], label='cos')\n", "ascii_art(cosxy)" ] }, { "cell_type": "code", "execution_count": 29, "id": "4f786401-ff63-48e6-ab14-df2d4c59366c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ " ___+____\n", " / / / \n", " cos sin -1\n", " | | \n", " *_ x \n", " / / \n", "x y " ] }, "execution_count": 29, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sinx = LabelledOrderedTree([Lx], label='sin')\n", "f2 = LabelledOrderedTree([cosxy, sinx, Lm1], label='+')\n", "ascii_art(f2)" ] }, { "cell_type": "code", "execution_count": 30, "id": "c60265b0-877c-46b4-b619-3403a8b9207d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ " _____*_______\n", " / / \n", " _^___ ___+____\n", " / / / / / \n", " _+__ -1 cos sin -1\n", " / / / | | \n", "x y -1 *_ x \n", " / / \n", " x y " ] }, "execution_count": 30, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tree = LabelledOrderedTree([f1, f2], label='*')\n", "ascii_art(tree)" ] }, { "cell_type": "markdown", "id": "7eed4cb1-2ae0-469e-bb55-95588e5b410f", "metadata": {}, "source": [ "# question 6" ] }, { "cell_type": "markdown", "id": "f7d3881a-e7a6-4835-a2c3-d14762722ac7", "metadata": {}, "source": [ "Let $\\displaystyle \\frac{1}{n} \\sum_{k=0}^{n-1} \n", "\\cos\\left(\\left( \\frac{k}{n} \\right)^2\\right)\n", "\\approx \\int_0^1 \\cos(x^2) dx$\n", "\n", "where ``s = lambda n: float(sum([cos((k/n)**2) for k in range(n)]))/n``\n", "\n", "is a function which defines the approximation for the integral.\n", "\n", "1. Time the execution of ``s`` for $n = 10000$.\n", " Explain why ``s`` is inefficient.\n", "\n", "2. Apply vectorization to improve the efficiency. Verify the correctness.\n", "\n", " Time the execution of the vectorized function for $n = 10000$,\n", " compare with timings of ``s``." ] }, { "cell_type": "markdown", "id": "dd853506-88d7-4dea-8fb5-bf509d440f31", "metadata": {}, "source": [ "## answer to question 6" ] }, { "cell_type": "code", "execution_count": 31, "id": "1944b91f-595f-4660-935e-46eeeccbbbe1", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 32, "id": "b23dbea5-b088-494a-9d02-28dd52bd72b3", "metadata": {}, "outputs": [], "source": [ "s = lambda n: float(sum([cos((k/n)**2) for k in range(n)]))/n" ] }, { "cell_type": "code", "execution_count": 33, "id": "f3b56eab-d533-4535-ad3a-e27cc6a34f08", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "5 loops, best of 3: 650 ms per loop" ] }, "execution_count": 33, "metadata": {}, "output_type": "execute_result" } ], "source": [ "timeit('s(10000)')" ] }, { "cell_type": "markdown", "id": "2280b531-86ea-4d3c-a94e-ab79018da04d", "metadata": {}, "source": [ "The Python function is inefficient because all computations are done in exact rational arithmetic and the evaluation to a double happens only at the very end." ] }, { "cell_type": "code", "execution_count": 34, "id": "ef4a9d3c-7161-43cc-8985-6890af24fda1", "metadata": {}, "outputs": [], "source": [ "def vs(n):\n", " \"\"\"\n", " Vectorized version of s.\n", " \"\"\"\n", " from numpy import arange, cos, sum\n", " x = arange(0,n)/n\n", " return sum(cos(x^2))/n" ] }, { "cell_type": "markdown", "id": "d857ab07-c8c1-4929-98c7-08134f653cca", "metadata": {}, "source": [ "We first verify the correctness." ] }, { "cell_type": "code", "execution_count": 35, "id": "1536c21a-9539-4776-97ac-1f65041b163c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.9045472213825253" ] }, "execution_count": 35, "metadata": {}, "output_type": "execute_result" } ], "source": [ "s1 = s(10000)\n", "s1" ] }, { "cell_type": "code", "execution_count": 36, "id": "60a5d52d-7583-45e3-8ba7-7671b982a3d0", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.904547221382527" ] }, "execution_count": 36, "metadata": {}, "output_type": "execute_result" } ], "source": [ "s2 = vs(10000)\n", "s2" ] }, { "cell_type": "code", "execution_count": 37, "id": "33fb27ba-f572-4664-a809-1b73b84b4a39", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.6653345369377348e-15" ] }, "execution_count": 37, "metadata": {}, "output_type": "execute_result" } ], "source": [ "abs(s1-s2)" ] }, { "cell_type": "markdown", "id": "75b59c2a-e146-436b-a609-16258438bd93", "metadata": {}, "source": [ "As the difference is close to the machine precision, the vectorized version is correct." ] }, { "cell_type": "code", "execution_count": 38, "id": "34ff4c2c-d6b4-44da-a524-00a3e18d58bf", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "625 loops, best of 3: 52.8 μs per loop" ] }, "execution_count": 38, "metadata": {}, "output_type": "execute_result" } ], "source": [ "timeit('vs(10000)')" ] }, { "cell_type": "markdown", "id": "f19d987a-f6d5-4120-93c4-72bc64a41421", "metadata": {}, "source": [ "We see a dramatic difference in execution times, going from hundreds of milliseconds to 50 microseconds." ] }, { "cell_type": "markdown", "id": "78001295-19f2-424f-b6a5-743f363f06bc", "metadata": {}, "source": [ "# question 7" ] }, { "cell_type": "markdown", "id": "c4ae36bc-cd39-4356-b4be-530b44260ea7", "metadata": {}, "source": [ "Let $p = 3 x^7 - 8 x^6 + 5 x^5 + 2 x^4 - 7 x^3 + 9 x^2 + 6 x + 1$\n", "be an expression in $x$.\n", "\n", "1. What is the best form to numerically evaluate $p$ in the fastest way?\n", "\n", "2. Illustrate your answer with timings of the evaluation of this form\n", " at a random number, comparing with the evaluation of $p$,\n", " given as a symbolic expression." ] }, { "cell_type": "markdown", "id": "2eaec5a5-40f9-4fff-82fa-b8d0b2f187dd", "metadata": {}, "source": [ "## answer to question 7" ] }, { "cell_type": "code", "execution_count": 39, "id": "b4d4a624-6d1b-4e32-a68a-0d2d3938148c", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "markdown", "id": "2bba86e3-2d8e-4e65-ad8a-2ccddc46266b", "metadata": {}, "source": [ "The fastest way to evaluate $p$ is to make a fast callable object of the Horner form of $p$." ] }, { "cell_type": "code", "execution_count": 40, "id": "6a33c245-b15d-411a-ae76-46fe5c99c90c", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle 3 \\, x^{7} - 8 \\, x^{6} + 5 \\, x^{5} + 2 \\, x^{4} - 7 \\, x^{3} + 9 \\, x^{2} + 6 \\, x + 1\\)" ], "text/latex": [ "$\\displaystyle 3 \\, x^{7} - 8 \\, x^{6} + 5 \\, x^{5} + 2 \\, x^{4} - 7 \\, x^{3} + 9 \\, x^{2} + 6 \\, x + 1$" ], "text/plain": [ "3*x^7 - 8*x^6 + 5*x^5 + 2*x^4 - 7*x^3 + 9*x^2 + 6*x + 1" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x = var('x')\n", "p = 3*x^7 - 8*x^6 + 5*x^5 + 2*x^4 - 7*x^3 + 9*x^2 + 6*x + 1\n", "show(p)" ] }, { "cell_type": "code", "execution_count": 41, "id": "7dbce662-6c85-4886-ad3c-481eabb6d5ad", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle {\\left({\\left({\\left({\\left({\\left({\\left(3 \\, x - 8\\right)} x + 5\\right)} x + 2\\right)} x - 7\\right)} x + 9\\right)} x + 6\\right)} x + 1\\)" ], "text/latex": [ "$\\displaystyle {\\left({\\left({\\left({\\left({\\left({\\left(3 \\, x - 8\\right)} x + 5\\right)} x + 2\\right)} x - 7\\right)} x + 9\\right)} x + 6\\right)} x + 1$" ], "text/plain": [ "((((((3*x - 8)*x + 5)*x + 2)*x - 7)*x + 9)*x + 6)*x + 1" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "h = p.horner(x)\n", "show(h)" ] }, { "cell_type": "code", "execution_count": 42, "id": "b97fa189-b827-44dd-9fb7-908b69b4c68b", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[('load_arg', 0),\n", " ('load_const', 3),\n", " 'mul',\n", " ('load_const', -8),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', 5),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', 2),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', -7),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', 9),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', 6),\n", " 'add',\n", " ('load_arg', 0),\n", " 'mul',\n", " ('load_const', 1),\n", " 'add',\n", " 'return']" ] }, "execution_count": 42, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f = fast_callable(h, vars=['x'])\n", "f.op_list()" ] }, { "cell_type": "markdown", "id": "48c417f7-bfd0-4e4e-a9a1-2daefbc7acb9", "metadata": {}, "source": [ "Let us compare the evaluate at a random real number." ] }, { "cell_type": "code", "execution_count": 43, "id": "b3c00e08-0f05-487c-ad25-07397461338c", "metadata": {}, "outputs": [], "source": [ "z = RR.random_element()" ] }, { "cell_type": "code", "execution_count": 44, "id": "c6f5d196-69b2-4c67-84df-315f0239ab5e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "-0.575080005661158" ] }, "execution_count": 44, "metadata": {}, "output_type": "execute_result" } ], "source": [ "p(x=z)" ] }, { "cell_type": "code", "execution_count": 45, "id": "ce87c7e3-7e55-4872-a617-9966b8b6a505", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "-0.575080005661155" ] }, "execution_count": 45, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f(z)" ] }, { "cell_type": "code", "execution_count": 46, "id": "bb0011b0-488a-4190-acee-f9da626e36f0", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "625 loops, best of 3: 38.6 μs per loop" ] }, "execution_count": 46, "metadata": {}, "output_type": "execute_result" } ], "source": [ "timeit('p(x=z)')" ] }, { "cell_type": "code", "execution_count": 47, "id": "12fb0651-1ec2-4ed4-9c30-02c1f2d095af", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "625 loops, best of 3: 5.25 μs per loop" ] }, "execution_count": 47, "metadata": {}, "output_type": "execute_result" } ], "source": [ "timeit('f(z)')" ] }, { "cell_type": "markdown", "id": "96a3ea83-393b-444d-9600-d0d06e3b4246", "metadata": {}, "source": [ "We see that the fast callable evaluates almost 7 times faster than the expression." ] }, { "cell_type": "markdown", "id": "85e4387d-2569-49f2-b686-832403d99a5b", "metadata": {}, "source": [ "# question 8" ] }, { "cell_type": "markdown", "id": "08b3a3fd-7a4d-4091-b9e7-9633972c2256", "metadata": {}, "source": [ "Transform $p = (u+1)^6 - v^4$ into $r = ((u+1)^3 + v^2)((u+1)^3 - v^2)$,\n", "without typing $r$." ] }, { "cell_type": "markdown", "id": "336dcefd-3e04-400e-a4db-2e22f8533fc0", "metadata": {}, "source": [ "## answer to question 8" ] }, { "cell_type": "code", "execution_count": 48, "id": "9424445e-989d-4e9f-9ee7-030315455cd9", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 49, "id": "7546b3a8-81aa-43e1-b63f-fcd77c5d7b1f", "metadata": {}, "outputs": [], "source": [ "u, v = var('u, v')" ] }, { "cell_type": "code", "execution_count": 50, "id": "c15eb9cc-6903-4efa-af41-d7efd4875fc2", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle {\\left(u + 1\\right)}^{6} - v^{4}\\)" ], "text/latex": [ "$\\displaystyle {\\left(u + 1\\right)}^{6} - v^{4}$" ], "text/plain": [ "(u + 1)^6 - v^4" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "p = (u + 1)^6 - v^4\n", "show(p)" ] }, { "cell_type": "markdown", "id": "c58ef9dc-ed27-44bb-aaeb-7bfa2c920867", "metadata": {}, "source": [ "To prevent $(u+1)^6$ from expanding, we replace $u+1$ by $w$ in $p$." ] }, { "cell_type": "code", "execution_count": 51, "id": "e61fa38a-cc40-4ee7-b054-fbd22622b8d4", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle w^{6} - v^{4}\\)" ], "text/latex": [ "$\\displaystyle w^{6} - v^{4}$" ], "text/plain": [ "w^6 - v^4" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "w = var('w')\n", "pw = p.subs({u+1: w})\n", "show(pw)" ] }, { "cell_type": "code", "execution_count": 52, "id": "df6f2c94-0dcc-4d4b-8ec0-1b8d51ecdce1", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle {\\left(w^{3} + v^{2}\\right)} {\\left(w^{3} - v^{2}\\right)}\\)" ], "text/latex": [ "$\\displaystyle {\\left(w^{3} + v^{2}\\right)} {\\left(w^{3} - v^{2}\\right)}$" ], "text/plain": [ "(w^3 + v^2)*(w^3 - v^2)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fpw = factor(pw)\n", "show(fpw)" ] }, { "cell_type": "code", "execution_count": 53, "id": "4e6dae6f-5e69-42bb-8689-40bf1083d313", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle {\\left({\\left(u + 1\\right)}^{3} + v^{2}\\right)} {\\left({\\left(u + 1\\right)}^{3} - v^{2}\\right)}\\)" ], "text/latex": [ "$\\displaystyle {\\left({\\left(u + 1\\right)}^{3} + v^{2}\\right)} {\\left({\\left(u + 1\\right)}^{3} - v^{2}\\right)}$" ], "text/plain": [ "((u + 1)^3 + v^2)*((u + 1)^3 - v^2)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "r = fpw(w = u+1)\n", "show(r)" ] }, { "cell_type": "markdown", "id": "38bce6f8-78e2-4154-bc1d-2447291631d9", "metadata": {}, "source": [ "# question 9" ] }, { "cell_type": "markdown", "id": "16bfb520-0c9e-4052-8e98-a824d946da14", "metadata": {}, "source": [ "Consider the expressions\n", "$\\displaystyle q = \\frac{x^2 + 4 x - 12}{x-2}$ and $p = x + 6$.\n", "\n", "1. Are the expressions *symbolically* the same?\n", "\n", "2. Are $p$ and $q$ *numerically* equivalent?\n", "\n", "Justify your answers by the appropriate computations." ] }, { "cell_type": "markdown", "id": "dd38888f-b774-4850-b45c-211c6d9969c4", "metadata": {}, "source": [ "## answer to question 9" ] }, { "cell_type": "code", "execution_count": 54, "id": "bc5b5a05-69a5-4cc2-bf84-f978cc05f56a", "metadata": {}, "outputs": [], "source": [ "reset()" ] }, { "cell_type": "code", "execution_count": 55, "id": "b07daede-6a47-4a36-86e8-9cb55fa64357", "metadata": {}, "outputs": [], "source": [ "x = var('x')" ] }, { "cell_type": "code", "execution_count": 56, "id": "55fa221e-c56b-4920-a30c-7daf28760df6", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\frac{x^{2} + 4 \\, x - 12}{x - 2}\\)" ], "text/latex": [ "$\\displaystyle \\frac{x^{2} + 4 \\, x - 12}{x - 2}$" ], "text/plain": [ "(x^2 + 4*x - 12)/(x - 2)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "q = (x^2 + 4*x - 12)/(x-2)\n", "show(q)" ] }, { "cell_type": "code", "execution_count": 57, "id": "3c7ce2fd-8cae-4b0a-9b16-90a07cb59663", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle x + 6\\)" ], "text/latex": [ "$\\displaystyle x + 6$" ], "text/plain": [ "x + 6" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "p = x + 6\n", "show(p)" ] }, { "cell_type": "code", "execution_count": 58, "id": "6baf13b2-4017-43fe-8f35-6d76bc40256f", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle x + 6\\)" ], "text/latex": [ "$\\displaystyle x + 6$" ], "text/plain": [ "x + 6" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "nq = q.normalize()\n", "show(nq)" ] }, { "cell_type": "code", "execution_count": 59, "id": "9c1b3e2d-4612-436d-a8b8-edacb0e3f818", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0" ] }, "execution_count": 59, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nq - p" ] }, { "cell_type": "markdown", "id": "0be1dd33-8b8a-4114-b8b5-dfea02e5cc5d", "metadata": {}, "source": [ "When $q$ is normalized, the common factor is removed and the normalized expression equals $p$. Symbolically, both expressions are thus they same." ] }, { "cell_type": "markdown", "id": "d343b2fb-f110-422c-adbc-c095856719dc", "metadata": {}, "source": [ "The numerical equality test compares the two expressions, evaluated at a random number." ] }, { "cell_type": "code", "execution_count": 60, "id": "d7d1ac22-5082-4d53-a72a-9d54589c5c14", "metadata": {}, "outputs": [], "source": [ "z = RR.random_element()" ] }, { "cell_type": "code", "execution_count": 61, "id": "4e03ff15-5ce8-4314-9dd8-4803c8f137ad", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "5.61039556909972" ] }, "execution_count": 61, "metadata": {}, "output_type": "execute_result" } ], "source": [ "qz = q(x=z)\n", "qz" ] }, { "cell_type": "code", "execution_count": 62, "id": "2faf671a-9f26-4cde-ae3f-1228ab595609", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "5.61039556909972" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pz = p(x=z)\n", "pz" ] }, { "cell_type": "code", "execution_count": 63, "id": "519f33f4-b9ba-4466-b399-bc3a3da168b8", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.000000000000000" ] }, "execution_count": 63, "metadata": {}, "output_type": "execute_result" } ], "source": [ "abs(qz-pz)" ] }, { "cell_type": "markdown", "id": "ae4f6f96-4847-4377-a27b-d6fc6385d67e", "metadata": {}, "source": [ "Numerically, we see that both expressions yield the same function values at random points. Therefore, also numerically $p$ and $q$ are the same." ] }, { "cell_type": "markdown", "id": "48344812-ad66-4bdb-9824-b51c42156edc", "metadata": {}, "source": [ "Although we know that ``q(x=2)`` leads to a division by zero, we can still say that the two expressions are identical wherever their domain is." ] } ], "metadata": { "kernelspec": { "display_name": "SageMath 10.3", "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.10.14" } }, "nbformat": 4, "nbformat_minor": 5 }