import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

def proche(y) :
    a = np.floor(y)
    b = a+1
    if y-a < 0.5 :
        return a
    else :
        return b

def fact(n) :
    r = 1
    for k in range(1,n+1) :
        r *= k
    return r
        
def a(n) :
    return proche(fact(n)*np.exp(-1))
    
def b(n) :
    s = 0
    for k in range(n+1) :
        s += (-1)**k/fact(k)
    return s * fact(n)
    
# for k in range(1,17) :
#     print(a(k),b(k))
    
def d(n) :
    return fact(n)*np.exp(-1)-b(n)
    
def J(n) :
    def fct(x) :
        return x**n*np.exp(x)
    return quad(fct,0,1)[0]
    
# for k in range(1,17) :
#     print(d(k),J(k)/d(k))
    
X,Y = [],[]
for k in range(1,17) :
    X.append(k)
    Y.append(d(k))
plt.plot(X,Y)