#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Wed May 15 21:20:31 2024

@author: matthieu
"""

#%%sujet 1
import numpy.linalg as la
import numpy as np
def K(n):
    M=np.zeros([n+1,n+1])
    for i in range(0,n): # attention décalage
        M[i][i+1]=i+1

    for j in range(0,n):
        M[j+1][j]=-n-1+j+1
    return M

for n in range(1,11):
    print(la.eigvals(K(n)))

#%%
## Sujet 2
import random as rd

def simul(n):
    return [rd.randint(1,n) for i in range(n)]

def X(n):
    return (len(set(simul(n))))

def moyenne(n,N=10**5):
    S=0
    for i in range(N):
        S+=X(n)
    return S/N
#%% Sujet 3
##
import numpy as np
import numpy.linalg as la
import matplotlib.pyplot as plt
A=np.array([[2,1,1],[1,2,1],[0,0,3]])
print(la.eig(A))
P=np.array([[1,0,1],[1,0,-1],[0,1,0]])
Pi=la.inv(P)           
a=0
b=1
Npas =100
x=np.linspace(a,b,Npas)
y=np.exp(3*x)
plt.plot(x,y,color='goldenrod')

##
#% sujet 4
import numpy.random as rd



def expo(lam):
    return -np.log(rd.random())/lam
def mini(lam,mu):
    return min(expo(lam),expo(mu))

def simulation(N):
    nbC2,nbC3=0,0
    for i in range(N):
        X1,X3=expo(1),expo(1)
        X2,X4=expo(2),expo(2)
        if X1>=X2 and X3>=X2:
            nbC2+=1
        if X2>=X3 and X4>=X3:
            nbC3+=1
    return nbC2/N,nbC3/N
        
##
#%% sujet 5
import random as rd
def binom(n,p):
    NbSucces=0
    for i in range(n):
        if rd.random()<p:
            NbSucces+=1
    return NbSucces

def Y(n):
    R=[]
    k=1
    for i in range(n):
        k=binom(3, k/3)
        R.append(k)
    return R
def EsperanceY(n,N):
    R=np.zeros(n)
    k=1
    for  e in range(N):
        k=1
        for i in range(n):
            k=binom(3, k/3)
            R[i]+=k
    return R/N

##
#%% Sjet 6
import numpy.random as rd
import numpy as np
from scipy.special import binom, factorial 
def achat(mu,lam):
    N=rd.poisson(mu)
    X=rd.poisson(lam,N)
    return sum(X)
    
def f(n,x):
    F=[np.exp(x)-1]
    for k in range (0,n):
        S=x
        for i in range(k+1):
            S+=x*binom(k,i)*F[i]
        F.append(S)
    return F[n]

def proba(n,lam,mu):
    if n==0:
        return np.exp(lam*(np.exp(-mu)-1))
    else:
        return np.power(mu,n)*np.exp(-lam)*f(n,lam*np.exp(-mu))/factorial(n)
    
    
##
#%% sujet 7
import numpy as np
import matplotlib.pyplot as plt 
def dicho(f,a,b,y,eps):
    """ on suppose que f(x)-y change de signe entre a  et b"""
    while b-a>eps:
        c=(a+b)/2
        if (f(a)-y)*(f(c)-y)<0:
            b=c
        else:
            a=c
    return a

def f(x):
    return np.log(x)+x

Ln=[]
Lxn=[]
for i in range(1,11):
   Ln.append(i)
   Lxn.append(dicho(f,1,i,i,10**-3)) 
plt.axis('equal')
plt.xlim(0,10)
plt.ylim(0,10)
plt.plot(Ln,Lxn,'*',color='teal')
plt.show()


##
#%% Sujet 8
import numpy.random as rd
def Nbmin(L):
    """ renvoie le nombre d'occurences de la valeur minimale
    """
    Vmin =L[0]  #valeur du min
    Omin=1      # nd de min
    for i in range(1,len(L)):
        if L [i]<Vmin:
            Vmin=L[i]
            Omin=1
        elif L[i]==Vmin:
            Omin+=1
    return Omin
 

def geom(n,p):
    """ simulation de n va suivant une loi géometrique"""
    R=[]
    for i in range(n):
        r=1
        while rd.random()>p:
            r+=1
        R.append(r)
    return R
    
    
def N(n,p):
    return Nbmin(geom(n,p))



def estimation(n,p,NE):
    S=0
    for i in range(NE):
        S+=N(n,p)
    return S/NE
#test
Ntest=10**5
ListeP=[k/10 for k in range(1,10)]
Listen=[2**n for n in range(5)]
for p in ListeP:
    for n in Listen:
        q=1-p
        print(estimation(n,p,Ntest)-n*p/(1-q**n))



##
#%% Sujet 9
import numpy as np
import numpy.linalg as al
A=np.array ([[5/2,1,1/2],[1,2,1],[1/2,1,5/2]])
I=np.eye (3)
r=al.matrix_rank(A-I)
s=al.matrix_rank(A-2*I)
t=al.matrix_rank(A-4*I)

al.eig(A)
##
#%% Sujet 10
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm 

def f(x): # fonction de survie
    return norm.sf(x)

def phi(x):
    return np.exp(-x**2/2)/np.sqrt(2*np.pi)

def M(x):
    return f(x)/phi(x)

def fonc1(x):
    return 1/x
    
def fonc2(x):
    return 1/x-1/(x**3)

a=0.05
b=5
Npas=10**5
pas =(b-a)/Npas
X=[a+k*pas for k in range(1,Npas+1)]
YM=[M(x) for x in X]
plt.plot(X,YM,label='y=f(x)',color='teal')
Y1=[fonc1(x) for x in X]
plt.plot(X,Y1,label='y=1/x',color='goldenrod')
Y2=[fonc2(x) for x in X]
plt.plot(X,Y2,label='y=1 /x-1/x**3',color='orchid')
plt.rcParams['figure.dpi'] = 500
plt.legend()
plt.axis('equal')
plt.ylim((0,1.5))
plt.savefig('Sujet10.png')
plt.show()


##
#%% Sujet 11
import numpy as np


def integre(n,N):
    pas =(1-0)/N
    I=0
    for k in range(N):
        I+=((1-k*pas)*np.exp(k*pas))**n
    return I/N


N=10**5

for i in range(10):
    n=2**i
    print(np.sqrt(n)*integre(n,N))

##
#%% Sujet 12


def matrice(L):
    try:
        [a1,a2,a3]=L
    except:
        print("la liste n'a pas le bon format")
        
    else:
        return ([[1,a1,a1**2],[1,a2,a2**2],[1,a3,a3**2]])


def quantites(i,L):
    i=i-1  #décalage python/maths
    S=0
    P=1
    D=1
    for j in range(len(L)):
        if j!= i:
            S+=L[j]
            P*=L[j]
            D*=(L[i]-L[j])
    return [S,P,D]            
    

def inverse(L):
    try:
        [a1,a2,a3]=L
        s1,p1,d1=quantites(1,L)
        s2,p2,d2=quantites(2,L)
        s3,p3,d3=quantites(3,L)
        
    except:
        print("la liste n'a pas le bon format")
        
    else:


        return ([[p1/d1,p2/d2,p3/d3],[-s1/d1,-s2/d2,-s3/d3],[1/d1,1/d2,1/d3]])
    

import numpy.linalg as la
print(la.inv(matrice([2,3,4])))
print(inverse([2,3,4]))
##
#%% Sujet 13
from random import random
from math import log
def expo():
	return -log(random())

for i in range(10):
    N=10**i
    S=0
    for i in range(N):
       S+=1/max(expo(),expo())
    print(S/N)
##
#%% Sujet 13
import random as rd
import math as m
def simule(n):
    S=0
    for i in range(n):
        S+= rd.gauss(0,1)**2
    return S

n=100
N=10**3 # nombre d'expéreience
for t in [1,2,3]:
    seuil= 2*m.sqrt(n*t)+2*t
    NbSucces=0
    for i in range(N):
        if simule(n)-n>= seuil:
            NbSucces+=1
    print("A t'on", NbSucces/N, ' <', m.exp(-t),'?')
            

##
#%% SUjet 15

import numpy as np
def test(u,w):
    [a,b,c]=u
    [x,y,z]=w
    if abs(x**2+y**2+z**2-1)<10**-5:
        return([-b*y+a*z,b*a-c*z,-a*x+c*b])
    else:
        return False
        
## 
#%% Sujet 16
def Ndiff(L):
    #avec dictionnaire
    #facile à modifier 
    Occ={}
    for x in L:
        if x in Occ:
            Occ[x]+=1
        else:
            Occ[x]=1
    return len(Occ)
    
import  random as rd
def Z(N,k):
    L=[]
    for i in range(k):
        L.append(rd.randint(1,N))
    return (Ndiff(L))

def Z_bis(N,k):
    Occ={}
    for i in range(k):
        r=rd.randint(1,N)
        if r in Occ:
            Occ[r]+=1
        else:
            Occ[r]=1
    return len(Occ)

def esperance(N,k,Nbexp=10**3):
    S=0
    for i in range(Nbexp):
        S+=Z_bis(N,k)
    return S/Nbexp
##
#%%
#cas 1
#cas 3
N=10
Lk=[10**i for i in range(1,5)]
for k in Lk:
    print(esperance(k,N))


#cas 2
k=10
LN=[10**i for i in range(1,5)]
for i in LN:
    print(esperance(k,N))
##
##%cas 3
L=[10**i for i in range(1,5)]
for i in LN:
    print(esperance(i,i)/i)

import numpy as np
print(1-np.exp(-1))