Informatique/SPE/IPT/TP0 METHODES NUMERIQUES/main.py

198 lines
7.4 KiB
Python
Raw Normal View History

2020-09-01 17:53:50 +02:00
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")
2020-09-11 18:31:23 +02:00
if False:
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))
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
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]))
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
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])
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
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])
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
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])
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
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
ne = 100
n = 10
plt.plot(np.linspace(0, 1, num=n+1), euler(expo_euler, 0, 1, 1, n), 'bo')
plt.plot(np.linspace(0, 1, num=ne+1), np.exp(np.linspace(0, 1, num=ne+1)), 'r')
plt.show()
2020-09-01 17:53:50 +02:00
2020-09-11 18:31:23 +02:00
print("--- Calcul Matriciel ---\n\n")
with Question(1):
with Question(2, 2):
with Question(1,3):
def matrice_nulle(n, m):
return [[0 for _ in range(m)] for _ in range(n)]
with Question(2,3):
def dimension(A):
if len(A):
return len(A),len(A[0])
else:
return 0
with Question(3,3):
def addition(A, B):
n,m = dimension(A)
return [[A[i][j] + B[i][j] for j in range(m)] for i in range(n)]
A3=[[1,2,3],
[0,1,0],
[-1,2,-6]]
B3=[[-1,5,0],
[1,0,-2],
[1,2,3]]
print(addition(A3, B3), '\n', np.array(A3)+np.array(B3))
print("Le résultat pour cette exemple est le même.")
with Question(4,3):
def transposee(A):
n,m = dimension(A)
return [[A[j][i]for j in range(n)] for i in range(m)]
print(transposee([[1,2,3],
[4,5,6]]))
print(list(map(list,zip(*[[1,2,3],[4,5,6]])))) # Parce que c'est amusant de chercher à faire une version courte.
with Question(5,3):
def multiple(A, coef):
n,m = dimension(A)
return [[coef * A[i][j] for j in range(m)] for i in range(n)]
print(multiple([[1,2],
[3,4]], 2))
with Question(6,3):
def multiplication(A, B):
n,p = dimension(A)
_,m = dimension(B)
return [[sum(A[i][k] * B[k][j] for k in range(p))for j in range(m)] for i in range(n)]
A6 = [[1,2,3],
[4,5,6]]
B6 = [[3,6],
[2,5],
[1,4]]
print(multiplication(A6,B6))
with Question(7,3):
def puissance(A, n):
_,m = dimension(A)
res = [[1 if j==i else 0 for j in range(m)] for i in range(m)]
for _ in range(0, n):
print(n)
res = multiplication(res, A)
return res
print(puissance([[0,1],
[0,0]], 2))
with Question(3, 2):
with Question(2,3):
with Question(1,4):
def echange_ligne(A, i, j):
A[i], A[j] = A[j], A[i]
A1 = [[1,2,3],
[4,5,6]]
echange_ligne(A1, 0, 1)
print(A1)
with Question(2,4):
def transvection(A, i, j, mu):
A[i] = addition([A[i]], multiple([A[j]], mu))[0]
transvection(A1, 0, 1, -1)
print(A1)
with Question(3, 4):
def pivot_partiel(A, j0):
n, m = dimension(A)
if n<j0 or m<j0:
raise IndexError("j0 trop grand")
maxi = A[j0][j0]
for i in range(j0, n):
maxi = max(A[i][j0], maxi)
return maxi
print(pivot_partiel([[1,2,3],
[4,5,6],
[7,3,9]], 1))
with Question(4,4):
def devient_triangle(A):
pass