diff --git a/0.3 Some Exercises.ipynb b/0.3 Some Exercises.ipynb
index 6747bf2..cd6d459 100644
--- a/0.3 Some Exercises.ipynb
+++ b/0.3 Some Exercises.ipynb
@@ -1,418 +1,418 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Some exercises\n",
"\n",
"**Course**: [_Analyse Numérique pour SV_](https://moodle.epfl.ch/course/info.php?id=) (MATH-2xx)\n",
"\n",
"Original by prof. Fabio Nobile\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"- Before starting to work on this exercise session, you should review the material in the short `Python tutorial` available on Moodle."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 0 -- Learning Python using Python."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Type the following commands and look at the results. Begin by importing numpy and matplotlib.pyplot as "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"These are the libraries that python relies upon for scientific computing and plotting. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"| No. | Code | Instruction |\n",
"| --- | --- | --- |\n",
"| 1 | `a=7` | Define $a$ as a scalar variable which has type `int` |\n",
"| 2 | `a=7.1` | Define $a$ as a scalar variable which has type `float` |\n",
"| 3 | `b=np.array([1,2,3])` | Define $b\\in\\mathbb{R}^{1 \\times 3}$. What happens if you replace `,` by `''` or `;`? |\n",
"| 4 | `g=b[2]` | Return the **third** element of $b$. Note: Indexation starts at 0. |\n",
"| 5 | `E=np.array([[1, 2, 3],[4, 5, 6]])` | Create matrix $E\\in\\mathbb{R}^{2\\times3}$. Note: Python is case sensitive. |\n",
"| 6 | `E[0,:]` | Return first row of $E;$ $(1,4)^T$ |\n",
"| 7 | `h=E[1,2]` | Return the element $E_{2,3}$ of the matrix $E$ ($E_{2,3}=6$). |\n",
"| 8 | `%whos` | Fetch information about the used variables. |\n",
"| 9 | `help(np.sin)` | Fetch help for the command ` np.sin`. Same as `? np.sin` in iPython. |\n",
"| 10 | `A=np.eye(3)` | Return the identity matrix $3\\times 3$ |\n",
"| 11 | `I=np.ones((4,4))` | Return matrix $4\\times 4$ for which each value is 1 |\n",
"| 12 | `F=np.zeros((2,3))` | Return $F\\in\\mathbb{R}^{2\\times 3}$ for which each value is 0 |\n",
"| 13 | `F.T` | Perform transpose of $F$ |\n",
"| 14 | `x=np.linspace(0,1,3)` | Return a vector of length $3$ whose values are equally distributed between $0$ and $1$ |\n",
"| 15 | `D=np.diag(b)` | Return diagonal matrix with diagonal given by vector b. |\n",
"| 16 | `np.shape(F)` | Return the number of lines and columns of $F$ in a line vector. |\n",
"| 17 | `A@D` | Return product of two matrices |\n",
"| 18 | `E@x` | Return matrix vector product |\n",
"| 19 | `A[-1,1]` | Return element of $A$ which is on the last line, second column |\n",
"| 20 | `H=np.array([[1, 3],[9, 11]]);np.linalg.solve(H,E)` | Compute $H^{-1}E$. Note: Never use `np.linalg.inv(H)@E` |\n",
"| 21 | `np.linalg.det(H)` | Compute determinant of $H$ |"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"## type here and execute the 21 commands above\"\n",
"#1\n",
"a=7\n",
"print('#1',a),\n",
"#2 continue by yourself"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 1\n",
"\n",
"Fibonacci's series is defined as follows:\n",
"$$\\begin{align*}\n",
"F_0&=0\\\\\n",
"F_1&=1\\\\\n",
"F_{n}&=F_{n-1}+F_{n-2}, \\quad \\forall n \\geq 2\n",
"\\end{align*}$$\n",
"* Complete the `Python` function below which computes the $30^{\\text{th}}$ Fibonacci's number.\n",
"\n",
"* Write a function which takes as an argument a natural number $n$ and returns the $n^{\\text{th}}$ Fibonacci's number.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"#COMPLETE THIS CODE\n",
"def Fibonacci(n):\n",
" # n should be a positive integer\n",
" # Complete this code\n",
" Fn=None #Change this value\n",
" # ...\n",
" return Fn\n",
"\n",
"F30=Fibonacci(30)\n",
"print(F30)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 2\n",
"Let us consider the matrices \n",
"\\begin{align*}\n",
"A=\n",
"\\begin{bmatrix}\n",
"5 & 3 & 0\\\\\n",
"1 & 1 & -4\\\\\n",
"3 & 0 & 0\n",
"\\end{bmatrix}, \\quad\n",
"B=\n",
"\\begin{bmatrix}\n",
"4 & 3 & 2\\\\\n",
"0 & 1 & 0\\\\\n",
"5 & 0 & 1/2\n",
"\\end{bmatrix}.\n",
"\\end{align*}\n",
"Compute (without using any loops), the matrix $C=AB$ (matrix product) and the matrix $D$ which\n",
"has as elements $D_{ij}=A_{ij}B_{ij}$ (element-wise product)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Solution:\n",
"\n",
"We can use the Python commands `@` for matrix multiplication and `*` for element-wise multiplication, as shown below:\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"A=np.array([ [5,3,0],\n",
" [1,1,-4],\n",
" [3,0,0]])\n",
"B=np.array([ [4,3,2],\n",
" [0,1,0],\n",
" [5,0,1/2]])\n",
"\n",
"\n",
"# Complete this code\n",
"C=None #Change this value\n",
"D=None #Change this value\n",
"\n",
"\n",
"print('C matrix')\n",
"print(C)\n",
"\n",
"print('D matrix')\n",
"print(D)\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 3\n",
"\n",
"Define (without using any loops) the\n",
"bidiagonal matrix of size $n=5$ whose main diagonal is a vector\n",
"of equally distributed points between $3$ and $6$, i.e. $$D=(3, 3.75, 4.5, 5.25, 6)^T,$$\n",
"\n",
"and the sub-diagonal is the vector of equally distributed points between $2$ and $3.5$, i.e. $$S=(2, 2.5, 3, 3.5)^T.$$\n",
"\n",
"\n",
"Tip: See ``https://numpy.org/doc/stable/reference/generated/numpy.diag.html``\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Solution\n",
"\n",
"Following the reference for the`diag` method of numpy we obtain:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Complete this code\n",
"D=None #Change this value\n",
"S=None #Change this value\n",
"M=None #Change this value\n",
"\n",
"print('The matrix M is')\n",
"print(M)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Remember that you can also obtain help from a python function using the command \n",
"`help(NAME OF FUNCTION)`. In this case: `help(np.diag)`"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"help(np.diag)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 4\n",
"Let us consider the vectors $x=(1,4,7,2,1,2)$ and $y=(0,9,1,4,3,0)$. \n",
"Compute (without using any loops for points a and b):\n",
"1. the product, component by component, between two vectors $x$ and $y$ (tip: use the operator `*`)\n",
"2. the scalar product between the same vectors $x$ and $y$ (tip: use the operator `@`)\n",
"3. a vector whose elements are defined by:\n",
"$v_1=x_1\\,y_n, \\quad v_2=x_2\\, y_{n-1},\\quad \\dots, \\quad v_{n-1}=x_{n-1}\\, y_2,\\quad v_n=x_n\\,y_1$.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Solution"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x=np.array([1,4,7,2,1,2])\n",
"y=np.array([0,9,1,4,3,0])\n",
"\n",
"# Complete this code\n",
"v1=None #Change this value\n",
"v2=None #Change this value\n",
"v3=None #Change this value\n",
"\n",
"print('v1 = '+str(v1))\n",
"print('v2 = '+str(v2))\n",
"print('v3 = '+str(v3))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 5\n",
"\n",
"Plot the functions\n",
"\\begin{align*}\n",
"f(x)=x^3 , \\quad x \\in [1,10]\n",
"\\end{align*}\n",
"\\begin{align*}\n",
"g(x)=\\exp(4x) , \\quad x \\in [1,10]\n",
"\\end{align*}\n",
"\n",
"in linear scale (i.e. using the function `plt.plot()` and logarithmic scale (i.e. using the function `plt.semilogy()` and `plt.loglog()`. We will use the definition interval to plot these graphs considering 200 equally distributed points.\n",
"\n",
"\n",
"\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Solution \n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Complete this code\n",
"f=None #Change this value to define inline function f\n",
"g=None #Change this value to define inline function g \n",
"\n",
"x=np.linspace(1,10,200)\n",
"\n",
"\n",
"#plot ins plot plot \n",
"plt.title(' Normal scale')\n",
"plt.plot(x,f(x),label=r'$f(x)$')\n",
"# Complete this code\n",
"# add a plot here for g\n",
"plt.legend()\n",
"plt.show()\n",
"\n",
"#plot in semilogy \n",
"\n",
"\n",
"# Complete this code\n",
"\n",
"\n",
"\n",
"\n",
"#plot in log-log\n",
"\n",
"# Complete this code\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercise 6\n",
"\n",
"1. Given \n",
"\\begin{align*}\n",
"f(x)=\\frac{x^2}{2}\\sin(x) , \\quad x \\in [1,20]\n",
"\\end{align*}\n",
"plot the function $f$ using 10, 20 and 100 equidistant points in the given interval.\n",
"Plot the three graphs on the same figure with three\n",
"different colors. Which one gives the best representation of $f$?\n",
"\n",
"2. Do the same for the functions:\n",
"\\begin{align*}\n",
"g(x) &= \\frac{x^3}{6}\\cos(\\sin(x))\\exp(-x)+\\left(\\frac{1}{1+x}\\right)^2 , \\quad x \\in [1,20] \\\\\n",
"h(x) &= x(1-x)+\\frac{\\sin(x)\\cos(x)}{x^3}, \\quad x \\in [1,20].\n",
"\\end{align*}\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### solution"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# complete the code here for question 1"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# complete the code here for question 2\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
- "version": "3.7.7"
+ "version": "3.11.5"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
diff --git a/1.1 Dichotomie.ipynb b/1.1 Dichotomie.ipynb
index 06fcc72..8968667 100644
--- a/1.1 Dichotomie.ipynb
+++ b/1.1 Dichotomie.ipynb
@@ -1,281 +1,281 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Equations non-linéaires\n",
"\n",
"**Objectif :** trouver les zéros (ou racines) d'une fonction $f:[a,b] \\rightarrow \\mathbf R$ : \n",
"$$\\alpha \\in [a,b] \\,:\\, f(\\alpha) = 0$$"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# importing libraries used in this book\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# defining the fonction\n",
"def f(x):\n",
" return x*np.sin(x*2.*np.pi) + 0.5*x - 0.25\n",
"\n",
"[a,b] = [-2,2]\n",
"\n",
"# points used to plot the graph \n",
"z = np.linspace(a, b, 100)\n",
"\n",
"plt.plot(z, f(z),'-')\n",
"\n",
"# labels, title, legend\n",
"plt.xlabel('x'); plt.ylabel('$f(x)$'); #plt.title('data')\n",
"plt.legend(['$f(t)$'])\n",
"plt.grid(True)\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1.1.1 Méthode de la bissection\n",
"\n",
""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Méthode de dichotomie ou bissection\n",
"\n",
"Si $f$ est continue et elle change de signe dans $[a,b]$, alors il existe au moins un $\\alpha$ tel que $f(\\alpha) = 0$.\n",
"\n",
"On peut alors définir l'algorithme suivant :\n",
"\n",
"$a^{(0)}=a$, $b^{(0)}=b$. Pour $k=0,1,...$\n",
"\n",
"1. $x^{(k)}=\\frac{a^{(k)}+b^{(k)}}{2}$\n",
"2. si $f(x^{(k)})=0$, alors $x^{(k)}$ est le zéro cherché.\n",
"\n",
" Autrement:\n",
"\n",
" 1. soit $f(x^{(k)})f(a^{(k)})<0$, alors le zéro\n",
" $\\alpha\\in[a^{(k)},x^{(k)}]$. \n",
"\n",
" On pose $a^{(k+1)}=a^{(k)}$ et $b^{(k+1)}=x^{(k)}$\n",
"\n",
" 2. soit $f(x^{(k)}f(b^{(k)})<0$, alors le zéro\n",
" $\\alpha\\in[x^{(k)},b^{(k)}]$. \n",
"\n",
" On pose $a^{(k+1)}=x^{(k)}$ et $b^{(k+1)}=b^{(k)}$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1.1.2 Bissection, critères d'arrêts\n",
"\n",
""
]
},
{
"cell_type": "markdown",
"metadata": {
"tags": []
},
"source": [
"## 1.1.3 Bissection, exemple\n",
"\n",
""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Exercice\n",
"Comprenez et completez la fonction suivante qui effectue l'algorithme de dichotomie \n",
"\n",
"Ensuite testez-la pour trouver la racine de la fonction $f(x) = x\\sin(2\\pi x) + \\frac12 x - \\frac14$ dans l'intervalle $[-1.5,1]$."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def bisection(a,b,fun,tolerance,maxIterations) :\n",
" # [a,b] interval of interest\n",
" # fun function\n",
" # tolerance desired accuracy\n",
" # maxIterations : maximum number of iteration\n",
" # returns:\n",
" # zero, residual, number of iterations\n",
" \n",
" if (a >= b) :\n",
" print(' b must be greater than a (b > a)')\n",
" return 0,0,0\n",
" \n",
" # what we consider as \"zero\"\n",
" eps = 1e-12\n",
" \n",
" # evaluate f at the endpoints\n",
" fa = fun(a)\n",
" fb = fun(b)\n",
"\n",
" if abs(fa) < eps : # a is the solution\n",
" zero = a\n",
" esterr = fa\n",
" k = 0\n",
" return zero, esterr, k\n",
" \n",
" if abs(fb) < eps : # b is the solution\n",
" zero = b\n",
" esterr = fb\n",
" k = 0\n",
" return zero, esterr, k\n",
"\n",
" if fa*fb > 0 :\n",
" print(' The sign of FUN at the extrema of the interval must be different')\n",
" return 0,0,0\n",
"\n",
"\n",
"\n",
" # We want the final error to be smaller than tol, \n",
" # i.e. k > log( (b-a)/tol ) / log(2) - 1\n",
"\n",
" nmax = int(np.ceil(np.log( (b-a)/tol ) / np.log(2))) - 1\n",
"\n",
" \n",
" # but nmax shall be smaller the the nmaximum iterations asked by the user\n",
" if ( maxIterations < nmax ) :\n",
" nmax = int(round(maxIterations))\n",
" print('Warning: nmax is smaller than the minimum number of iterations necessary to reach the tolerance wished');\n",
"\n",
" # vector of intermadiate approximations etc\n",
" x = np.zeros(nmax+1)\n",
"\n",
" # initial error is the length of the interval.\n",
" esterr = (b - a)\n",
"\n",
" # do not need to store all the a^k and b^k, so I call them with a new variable name:\n",
" ak = a\n",
" bk = b\n",
" # the values of f at those points are fa and fk\n",
"\n",
" for k in range(nmax+1) :\n",
"\n",
" # approximate solution is midpoint of current interval\n",
" # COMPLETE the code below\n",
" #> x[k] = \n",
" #> fx = \n",
" \n",
" # error estimator is the half of the previous error\n",
" #> esterr =\n",
"\n",
" # COMPLETE the code below\n",
"\n",
" # if we found the solution, stop the algorithm\n",
" if np.abs(fx) < eps :\n",
" # error is zero\n",
" #> zero = \n",
" #> esterr = \n",
" #> return \n",
"\n",
" if fx*fa < 0 : # alpha is in (a,x)\n",
" #> bk = \n",
" elif fx*fb < 0 : # alpha is in (x,b)\n",
" #> ak = \n",
" else :\n",
" error('Algorithm not operating correctly')\n",
"\n",
"\n",
"\n",
" zero = x[k];\n",
"\n",
" if esterr > tol :\n",
" print('Warning: bisection stopped without converging to the desired tolerance because the maximum number of iterations was reached');\n",
"\n",
" return zero, esterr, k, x \n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"[a,b] = [-1.5,1]\n",
"tol = 1e-10\n",
"maxIter = 4\n",
"zero, esterr, k, x = bisection(a,b,f,tol,maxIter)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from NonLinearEquationsLib import plotBisectionIterations"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plt.rcParams.update({'font.size': 16})\n",
"plt.figure(figsize=(8, 5))\n",
"\n",
"plotBisectionIterations(a,b,f,x)\n",
"\n",
"# plt.savefig('Bisection-iterations.png', dpi=600)\n",
"\n",
"print(f'The estimated root is {zero:2.12f}, the estimated error {esterr:1.3e} and the residual is {f(zero):1.3e}, after {k} iterations')"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
- "version": "3.7.7"
+ "version": "3.11.5"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
diff --git a/1.2 NetwonRaphson.ipynb b/1.2 NetwonRaphson.ipynb
index d7a0cc9..85f3ba2 100644
--- a/1.2 NetwonRaphson.ipynb
+++ b/1.2 NetwonRaphson.ipynb
@@ -1,271 +1,271 @@
{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# importing libraries used in this book\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1.2.1 Méthode de Newton\n",
"\n",
""
]
},
{
"cell_type": "markdown",
"metadata": {
"tags": []
},
"source": [
"## 1.2.2 Méthode de Newton, Graphique\n",
"\n",
""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Méthode de Newton\n",
"\n",
"Soit $f:\\mathbb R \\rightarrow\\mathbb R$ une fonction différentiable.\n",
"\n",
"Soit $x^{(0)}$ un point donné. On considère l'équation de la droite $y(x)$ qui\n",
"passe par le point $(x^{(k)},f(x^{(k)}))$ et qui a comme pente\n",
"$f'(x^{(k)})$,\n",
"\\begin{equation*}\n",
" y(x)=f'(x^{(k)})(x-x^{(k)})+f(x^{(k)}).\n",
"\\end{equation*}\n",
"On définit $x^{(k+1)}$ comme étant le point où cette droite\n",
"intersecte l'axe $x$, c'est-à-dire $y(x^{(k+1)})=0$. On en\n",
"déduit que :\n",
"$$ x^{(k+1)}=x^{(k)}-\\frac{f(x^{(k)})}{f'(x^{(k)})},\\,\\,\\,k=0,1,2\\ldots .$$\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Exercice (5, Série 1)\n",
"\n",
"On cherche les zéros de la fonction\n",
"\n",
"$$f(x) = \\dfrac{1}{2} \\sin \\left( \\dfrac{\\pi x}{2} \\right) + 1 - x \\; .$$\n",
"\n",
"\n",
"1. Vérifiez qu'il y a au moins un zéro $\\alpha$ dans l'intervalle $[0,2]$.\n",
"\n",
"2. Ecrivez la méthode de Newton pour trouver le zéro $\\alpha$ de la fonction $f(x)$ et\n",
"calculez la première itération à partir de la valeur initiale $x^{(0)}=1$.\n",
"\n",
"3. Calculez les zéros $\\alpha$ de la fonction\n",
" $f$ avec la méthode de Newton (fonction `newton` que vous devrez écrire)\n",
"\n",
"```python\n",
"def Newton( F, dF, x0, tol, nmax ) :\n",
" # NEWTON Find the zeros of a nonlinear equations.\n",
" # NEWTON(F,DF,X0,TOL,NMAX) tries to find the zero X of the \n",
" # continuous and differentiable function F nearest to X0 using \n",
" # the Newton method. DF is a function which take X and return the derivative of F.\n",
" # If the search fails an error message is displayed.\n",
" # \n",
" # returns the value of the\n",
" # residual R in X,the number of iterations N required for computing X and\n",
" # INC the increments computed by Newton.\n",
" \n",
" return x, r, n, inc\n",
"```\n",
"\n",
"Choisissez $x^{(0)} = 1$ comme point de départ pour la méthode\n",
"et utilisez une tolérance $tol=10^{-4}$ sur la valeur absolue de l'incrément\n",
"entre deux itérations\n",
"successives $|x^{(k+1)}-x^{(k)}|$. \n",
"\n",
"*Dans le cas de la méthode de Newton,\n",
"l'incrément est une bonne approximation de l'erreur.*\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def Newton( F, dF, x0, tol, nmax ) :\n",
" '''\n",
" NEWTON Find the zeros of a nonlinear equations.\n",
" NEWTON(F,DF,X0,TOL,NMAX) tries to find the zero X of the \n",
" continuous and differentiable function F nearest to X0 using \n",
" the Newton method. DF is a function which take X and return the derivative of F.\n",
" If the search fails an error message is displayed.\n",
" \n",
" Outputs : [x, r, n, inc, x_sequence]\n",
" x : the approximated root of the function\n",
" r : the absolute value of the residualin X\n",
" n : the number of iterations N required for computing X and\n",
" inc : the increments computed by Newton.\n",
" x_sequence : the sequence computed by Newton\n",
" '''\n",
" \n",
" # Initial values\n",
" n = 0\n",
" xk = x0\n",
"\n",
" # initialisation of loop components\n",
" # increments (in abs value) at each iteration\n",
" inc = []\n",
" # in case we wish to plot the sequence \n",
" x = [x0]\n",
" \n",
" # (Complete the code below)\n",
" \n",
" # Warning if not converged\n",
" if n > nmax :\n",
" print('Newton stopped without converging to the desired tolerance ')\n",
" print('because the maximum number of iterations was reached')\n",
" \n",
" # (Complete the code below)\n",
" \n",
" \n",
" return xk1, rk1, n, inc, np.array(x)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"f = lambda x : 0.5*np.sin(np.pi*x/2)+1-x\n",
"df = lambda x : 0.25*np.pi*np.cos(np.pi*x/2)-1\n",
"\n",
"x0 = 1\n",
"tol = 1e-4\n",
"nmax = 10\n",
"\n",
"zero, residual, niter, inc, x = Newton(f, df, x0, tol, nmax)\n",
"\n",
"print(f'The zero computed is {zero:1.4f}')\n",
"print(f'Newton stoppedconverged in {niter} iterations'); \n",
"print(f'with a residual of {residual:1.4e}.\\n');\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from NonLinearEquationsLib import plotNewtonIterations"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"[a,b] = [0,3.5]\n",
"x0 = 3\n",
"zero, residual, niter, inc, x = Newton(f, df, x0, tol, nmax)\n",
"\n",
"# Rezise plots, which are usually too small\n",
"plt.figure(figsize=(12, 4))\n",
"plt.rcParams.update({'font.size': 12})\n",
"\n",
"# Subplot 1 over 2, 1st one\n",
"plt.subplot(121)\n",
"\n",
"#plt.plot(range(MaxIterations), RelativeError, 'b:.')\n",
"plt.plot(range(niter), inc, 'b:.')\n",
" \n",
"plt.xlabel('n'); plt.ylabel('$\\\\delta x$');\n",
"plt.grid(True)\n",
"#plt.xscale('log') \n",
"plt.yscale('log')\n",
"plt.legend(['$|\\\\delta x|$'])\n",
"\n",
"\n",
"# Subplot 1 over 2, 2nd one\n",
"plt.subplot(122)\n",
"\n",
"plotNewtonIterations (a,b,f,x,200)\n",
"\n",
"plt.show()\n",
"\n",
"# Rezise plots, which are usually too small\n",
"plt.figure(figsize=(8, 4))\n",
"plt.rcParams.update({'font.size': 16})\n",
"\n",
"plotNewtonIterations (a,b,f,x,200)\n",
"# plt.savefig('Newton-iterations.png', dpi=600)\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"[a,b] = [-1,3.5]\n",
"z = np.linspace(a,b,200)\n",
"plt.plot(z,f(z), 'b-', x[6],f(x[6]), 'rx')\n",
"plt.ylabel('$f(x)$'); plt.xlabel('$x$');\n",
"\n",
"# Plot the x,y-axis \n",
"plt.plot([a,b], [0,0], 'k-',linewidth=0.1)\n",
"plt.plot([0,0], [np.min(f(z)),np.max(f(z))], 'k-',linewidth=0.1)\n",
"\n",
"plt.legend(['$f$','$\\\\alpha$'])\n",
"# plt.savefig('Newton-fx-alpha.png', dpi=600)\n",
"\n",
"plt.show()\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
- "version": "3.7.7"
+ "version": "3.11.5"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
diff --git a/2.2 Interpolation de Lagrange.ipynb b/2.2 Interpolation de Lagrange.ipynb
index 4faa0ff..acb3117 100644
--- a/2.2 Interpolation de Lagrange.ipynb
+++ b/2.2 Interpolation de Lagrange.ipynb
@@ -1,200 +1,200 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2.2 Interpolation de Lagrange\n",
"\n",
"### Base de Lagrange\n",
"\n",
"On considère les polynômes $\\varphi_k$, $k=0,\\dots, n$ de\n",
"degré $n$ tels que\n",
"\n",
"$$\\varphi_k(t_j)=\\delta_{jk}, \\qquad k,j=0,\\dots, n,$$\n",
"\n",
"où $\\delta_{jk}=1$ si $j=k$ et $\\delta_{jk}=0$ si $j\\neq k$.\n",
"Explicitement, on a\n",
"\n",
"$$\\varphi_k(t)=\\prod_{j=0,j\\ne k}^n\\dfrac{(t-t_j)}{(t_k-t_j)}.$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2.2.1 Base de Lagrange\n",
""
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# importing libraries used in this book\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice (théorique)** Vérifiez que\n",
"\n",
"1. $B = \\{\\varphi_k, k=0,\\dots, n\\}$ est une base de $P_n(\\mathbb R)$\n",
"2. Chaque polynôme $\\varphi_k$ est de degré $n$\n",
"3. $\\varphi_k(t_j)=\\delta_{jk}, \\qquad k,j=0,\\dots, n$\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Defining the Lagrange basis functions\n",
"def phi(t,k,z):\n",
" # the input variables are:\n",
" # t : the interpolatory points\n",
" # k : which basis function\n",
" # z : where to evaluate the function\n",
" \n",
"\n",
" # careful, there are n+1 interpolation points!\n",
" n = t.size - 1 \n",
"\n",
" # init result to one, of same type and size as z\n",
" result = np.zeros_like(z) + 1\n",
"\n",
" # first few checks on k:\n",
" if (type(k) != int) or (t.size < 1) or (k > n) or (k < 0):\n",
" raise ValueError('Lagrange basis needs a positive integer k, smaller than the size of t')\n",
" \n",
" # loop on n to compute the product\n",
" for j in range(0,n+1) :\n",
" if (j == k) :\n",
" continue\n",
" if (t[k] == t[j]) :\n",
" raise ValueError('Lagrange basis: all the interpolation points need to be distinct')\n",
"\n",
"# (Complete the code below)\n",
" result = \n",
"\n",
" return result\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Exemple\n",
"\n",
"Pour $n=2$, $t_0=-1$, $t_1=0$, $t_2=1$, les polynômes de la base de Lagrange\n",
"sont\n",
"\n",
"\\begin{equation*}\n",
"\\begin{array}{lcl}\n",
"\\varphi_0(t)&=&\\displaystyle\\dfrac{(t-t_1)(t-t_2)}{(t_0-t_1)(t_0-t_2)}=\n",
"\\dfrac{1}{2}t(t-1),\\\\[2mm]\n",
"\\varphi_1(t)&=&\\displaystyle\\dfrac{(t-t_0)(t-t_2)}{(t_1-t_0)(t_1-t_2)}=\n",
"-(t+1)(t-1),\\\\[2mm]\n",
"\\varphi_2(t)&=&\\displaystyle\\dfrac{(t-t_0)(t-t_1)}{(t_2-t_0)(t_2-t_1)}=\n",
"\\dfrac{1}{2}t(t+1).\n",
"\\end{array}\n",
"\\end{equation*}\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# plot the Lagrange Basis functions \n",
"t = np.linspace(-1., 1, 3)\n",
"z = np.linspace(-1.1, 1.1, 100)\n",
"\n",
"plt.plot(z, phi(t,0,z), 'g', z, phi(t,1,z), 'r', z, phi(t,2,z),':')\n",
"\n",
"plt.xlabel('t'); plt.ylabel('$\\\\varphi_{k}(t)$'); plt.title('Lagrange basis functions')\n",
"plt.legend(['$\\\\varphi_{0}$','$\\\\varphi_{1}$','$\\\\varphi_{2}$'])\n",
"plt.grid(True)\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Exercice\n",
"\n",
"Visualisez la base de Lagrange associée aux points $t_k = 1, 1.5, 2, 2.5, 3$.\n",
"Evaluez sur le graphique les valeurs de $\\varphi_k(t_j)$\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# plot the Lagrange Basis functions \n",
"t = np.linspace(1., 3, 5)\n",
"z = np.linspace(0.9, 3.1, 100)\n",
"\n",
"# (Complete the code below)\n",
"\n",
"\n",
"\n",
"\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Polynôme d'interpolation\n",
"\n",
"Le polynôme d'interpolation $\\Pi_n$ des\n",
"valeurs $y_j$ aux noeuds $t_j$, $j=0,\\dots,n$, s'écrit\n",
"\\begin{equation}\\label{e:int_lagr}\n",
" \\Pi_n(t)=\\sum_{k=0}^n y_k \\varphi_k(t),\n",
"\\end{equation}\n",
"car il vérifie $\\Pi_n(t_j)=\\sum_{k=0}^n y_k \\varphi_k(t_j)=y_j$.\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
- "version": "3.7.16"
+ "version": "3.11.5"
}
},
"nbformat": 4,
"nbformat_minor": 4
}