diff --git a/SPE/IPT/TP0 METHODES NUMERIQUES/main.py b/SPE/IPT/TP0 METHODES NUMERIQUES/main.py new file mode 100644 index 0000000..59ac9bc --- /dev/null +++ b/SPE/IPT/TP0 METHODES NUMERIQUES/main.py @@ -0,0 +1,107 @@ +import numpy as np + +class Question: + def __init__(self,i=1, l=0, n=''): + self.name=n + self.level=l + self.number=i + def __enter__(self): + print('\n' + (self.level*2)*' ' + f"-> {self.number}. : {self.name} -- Début") + return self + def __exit__(self, exc_type, exc_value, exc_traceback): + print((self.level*2)*' ' + f"<- {self.number}. : {self.name} -- Fin\n") + + +with Question(1): + with Question(1,1): + def dichotomie(f, n, a=0, b=2): + iterations = 0 + while abs(b-a) >= 2*10**(-n): + iterations += 1 + c = (a+b) / 2 + if f(c) == 0: + return c + elif f(a)*f(c) > 0: + a = c + else: + b = c + return c, iterations + f = lambda x: x**3 - 3*x**2 + 1 + print(dichotomie(f,5,0,1.5)) + print(dichotomie(np.sin, 3, 3, 4)) + print(dichotomie(np.sin, 10, 3, 4)) + + + with Question(2,1): + def newton(f, fp, x0, nbiter, epsilon=1e-10): + for i in range(nbiter): + x1 = x0 - f(x0) / fp(x0) + if abs(x1 - x0) <= epsilon: + break + x0 = x1 + return x1 + print(newton(lambda x: x**3 - 3*x**2 + 1, lambda x: 3*x**2 - 6*x, 1.5, 200, epsilon=1e-100)) + print(newton(np.sin, np.cos, 3, 100)) + print(abs(newton(np.sin, np.cos, 3, 100) - dichotomie(np.sin, 10, 3, 4)[0])) + + +with Question(2): + from math import cos, exp, sin, log + from scipy.integrate import quad + from time import time + with Question(1,2): + def rectangles(f, a, b, n): + res = 0 + for k in range(n): + res += f(a + k*(b - a)/n) + return res*(b - a)/n + def f1(x): + return x + def f2(x): + return x**10 + f3 = cos + f4 = exp + print(rectangles(f1, 0, 1, 1000), quad(f1, 0, 1)[0]) + print(rectangles(f2, 0, 1, 1000), quad(f2, 0, 1)[0]) + print(rectangles(f3, 0, np.pi/2, 1000), quad(f3, 0, np.pi/2)[0]) + print(rectangles(f4, -3, 3, 1000), quad(f4, -3, 3)[0]) + + + with Question(2,2): + def trapezes(f, a, b, n): + res = 0 + for k in range(n): + res += f(a + k*(b - a)/n) + f(a + (k+1)*(b - a)/n) + return res*(b - a)/(2*n) + print(trapezes(f1, 0, 1, 1000), quad(f1, 0, 1)[0]) + print(trapezes(f2, 0, 1, 1000), quad(f2, 0, 1)[0]) + print(trapezes(f3, 0, np.pi/2, 1000), quad(f3, 0, np.pi/2)[0]) + print(trapezes(f4, -3, 3, 1000), quad(f4, -3, 3)[0]) + + + with Question(3, 2): + def simpson(f, a, b, n): + res = 0 + for k in range(n): + res += f(a + k*(b - a)/n) + 4*f(((a+k*(b - a)/n) + (a + (k+1)*(b - a)/n))/2) + f(a + (k+1)*(b - a)/n) + return res*(b - a)/(6*n) + print(simpson(f1, 0, 1, 1000), quad(f1, 0, 1)[0]) + print(simpson(f2, 0, 1, 1000), quad(f2, 0, 1)[0]) + print(simpson(f3, 0, np.pi/2, 1000), quad(f3, 0, np.pi/2)[0]) + print(simpson(f4, -3, 3, 1000), quad(f4, -3, 3)[0]) + + +with Question(3): + import matplotlib.pyplot as plt + with Question(1, 2): + def euler(F, t0, y0, t1, n): + res = [y0]*(n+1) + for i in range(1,n+1): + res[i] = res[i-1] + (t1-t0)/n * F(t0 + i*(t1-t0)/n, res[i-1]) + return res + def expo_euler(t, yt): + return yt + n = 100 + plt.plot(np.linspace(0, 1, num=n+1), euler(expo_euler, 0, 1, 1, n)) + plt.plot(np.linspace(0, 1, num=n+1), np.exp(np.linspace(0, 1, num=n+1))) + plt.show()