#########################################
######## TP informatique bio-sép ########
########        2025-2026        ########
#########################################


#########################################
######### Thème 4 : Intégration #########
#########################################

import numpy as np
import matplotlib.pyplot as plt

## Q1 : schéma et explications
# La méthode des 'rectangles à gauche' approxime
# l'intégrale de f entre a et b par la somme des
# aires des rectangles construits sur les valeurs
# f(xk) donc la valeur de f 'à gauche' de Ik.
# La méthode des 'rectangles à droite' utilise
# f(xk+1) donc la valeur de f 'à droite' de Ik.

## Q2 : fonction 'Riemann' par la méthode
##      des rectangles à gauche

# voir script dans le poly de cours
def Riemann(f,a,b,n) :
    S = 0
    pas = (b-a)/n
    for k in range(n) :
        xk = a + k*pas
        S += f(xk)
    return S*pas


## Q3 : Test avec l'intégrale de la fonction identité

def id(x) :
    return x

for k in range(1,7) :
    print(Riemann(id,0,1,10**k))

print('\n')
print(Riemann(id,1,0,10**6))

# la fonction Riemann semble fonctionner aussi
# lorsque a > b.

## Q4 : Test avec x |--> 4 / (1+x²)

def g(x) :
    return 4/(1+x**2)

for k in range(1,7) :
    print(Riemann(g,0,1,10**k))

print('\n')
print(np.pi)

# On remarque que la convergence est assez lente :
# il faut 10**k calculs pour avoir une précision
# à 10**(-k) près.

## Q5 : fonction 'Riemann' par la méthode
##      des rectangles à droite

def RD(f,a,b,n) :
    S = 0
    pas = (b-a)/n
    for k in range(1,n+1) :
        xk = a + k*pas
        S += f(xk)
    return S*pas

## Q5bis : Comparaison RG et RD (tests avec g)

for k in range(1,7) :
    print('n = 10**(-', k, ')')
    print('   RG -> ',Riemann(g,0,1,10**k))
    print('   RD -> ',RD(g,0,1,10**k))

# sur cet exemple, RG donne une approximation par excès
# et RD une approximation par défaut.
# Les précisions sont comparables pour les 2 méthodes.

## Q6 : Programmation avec un pas passé en argument

def RG(f,a,b,pas) :
    S = 0
    x = a
    while x < b :
        S += f(x)
        x += pas
    return S*pas

## Q7 : ... et méthode des rectangles à droite

def RDpas(f,a,b,pas) :
    S = 0
    x = a + pas
    while x < b :
        S += f(x)
        x += pas
    return S*pas

## Q8 : Courbe de F définie par une intégrale

def F(x,n) :
    '''renvoie une valeur approchee de F(x)
    en utilisant n rectangles'''
    def f(t) :
        return 1/np.sqrt(t**4+t**2+1)
    return Riemann(f,x,2*x,n)

X = np.linspace(-10,10,400)
n = 10**4 # gare au temps de calcul si n est trop grand
Y = [F(x,n) for x in X]
plt.plot(X,Y)
plt.plot([-10,10],[0,0],color = 'black')
plt.plot([0,0],[-.4,.4],color = 'black')
plt.show()


## Q9 : Autre exemple

def G(x,n) :
    if x == 0 :
        return 0
    def g(t) :
        return np.sin(t*x)/(x+t)
    return Riemann(g,0,1,n)

X = np.linspace(0,30,400)
n = 10**4
Y = [G(x,n) for x in X]
plt.plot(X,Y)
plt.plot([0,30],[0,0],color = 'black')
plt.plot([0,0],[0,.3],color = 'black')
plt.show()