import matplotlib.pyplot as plt
import matplotlib as mpl
import numpy as np
from turtle import *
from random import random
import sys
from math import *
## Fonctions pour les piles
def estvide(p):
    return(p==[])

def newpile():
    return([])

def depile(p):
    return(p.pop())

def empile(x,p):
    p.append(x)

def sommet(p):
    return(p[-1])
##
##penup()
##speed(0)
##left(90)
##backward(300)
##pendown()
##
##def arbrec(n,l,th,c=0):
##    if n == 0 :
##        forward(l/3)
##        backward(l/3)
##    else :
##        tr = (random() - .5)/2 #Angles modifies entre -25 % et 25 % 
##        tr1 = (random() - .5)/2
##        lr = (random()-.25)*2/3 # Longueurs modifiees entre -17% et + 50%
##        lr1 = (random()-0.25)*2/3
##        pensize(5*n/(c+n))
##        tup = (0.5,0.3+0.7*(c-n)/(c+n),0)
##        pencolor(tup)
##        forward(l/3)
##        right(th)
##        arbrec(n-1,2*l*(1+lr)/3,th*(1+tr),c+1)
##        left(2*th)
##        arbrec(n-1,2*l*(1+lr1)/3,th*(1+tr1),c+1)
##        right(th)
##        penup()
##        backward(l/3)
##        pendown()
##
###arbrec(8,300,25,5)
##
##
##
##def trace(n):
##    for i in range(n):
##        forward(100/n**0.5)
##        left(360/n)
##
##
##def arbre(n,l,th):
##    if n == 0 :
##        forward(l)
##        backward(l)
##    else :
##        forward(l/3)
##        right(th)
##        arbre(n-1,2*l/3,th)
##        left(2*th)
##        arbre(n-1,2*l/3,th)
##        right(th)
##        backward(l/3)
##
##
##
##
##
##def Fc(a,b,L):
##    r = a%b
##    an = a//b
##    if r == 0:
##        L.append(an)
##        return(L)
##    else :
##        L.append(an)
##        return(Fc(b,r,L))
##
##
##
##def evalueFC(L):
##    if len(L) ==1 :
##        return(L[0])
##    else :
##        return(L[0]+1/evalueFC(L[1:]))
##
##def FcR(x,L,n):
##    an = int(x)
##    xn = x - an
##    if len(L) == n or abs(xn) < 10**(-10):
##        L.append(an)
##        return(L)
##    else :
##        L.append(an)
##        return(FcR(1/xn,L,n))
##
##plt.figure('Fraction')
##def Fce(x,n):
##    K = FcR(x,[],n)
##    T = [evalueFC(K[:i]) for i in range(1,len(K))]
##    plt.plot([x for i in T])
##    plt.plot(T)
##    plt.show()
##Fce(sqrt(exp(1)),8)

plt.figure('Percolation')
def probapond(a,b,p):
    """retourne a avec un proba de p et b avec une proba de 1-p, p dans [0,1]"""
    k = random()
    return(a*(k<p) + b*(k>p))


def cree_grille(n,p):
    sys.setrecursionlimit(900)
    M = np.zeros((n,n))
    for i in range(n):
        for j in range(n):
            M[i][j] = probapond(2,0,p)
    plt.imshow(M,interpolation='none', cmap='gist_stern', norm = mpl.colors.Normalize(vmin=0.,vmax=2.))
    plt.show()
    return(M)

M = cree_grille(500,0.6)
##def etat_case(M,i,j):
##    n = len(M)
##    if i >= 0 and j >= 0 and i < n and j < n:
##        return M[i][j]
##    else :
##        return 2


##def diffuse(M,a,b):
##    n = len(M)
##    L = [[-1,0],[0,-1],[1,0],[0,1]]
##    if etat_case(M,a,b) == 0 :
##        M[a][b] = 1
##        for i in L:
##            if etat_case(M,a+i[0],b + i[1]) ==0:
##                diffuse(M,a+i[0],b + i[1])
##        return M
##    else :
##        return M
##
##
##def remplit(M):
##    n = len(M)
##    for i in range(n):
##        diffuse(M,0,i)



def remplit2(M):
    n = len(M)
    p = newpile()
    for i in range(n):
        if M[0,i]==2:
            empile((0,i),p)
            M[0,i]=1
    while not estvide(p):
        a,b = depile(p)
        h = (a-1,b)
        d = (a,b+1)
        cb = (a+1,b)
        g = (a,b-1)
        if a >0 and M[h] == 2:
            M[h] = 1
            empile(h,p)
        if b < n-1 and M[d] == 2:
            M[d] = 1
            empile(d,p)
        if  a < n-1 and M[cb] == 2:
            M[cb] = 1
            empile(cb,p)
        if b>0 and M[g] == 2:
            M[g] = 1
            empile(g,p)
        
        
    
remplit2(M)
plt.imshow(M,interpolation='none', cmap='gist_stern', norm = mpl.colors.Normalize(vmin=0.,vmax=2.))
plt.show()


def percolation(p,n):
    M = cree_grille(n,p)
    remplit(M)
    i = 0
    while i < n and M[n-1][i] != 1: 
        i = i + 1
    return(i<n)



def Proba(p,n):
    P = 0
    for i in range(20):
        P = P + percolation(p,n)/20
    return P



##Lp = np.linspace(0,1,40)
##Ln = [i for i in range(1,40)]
##Y = [Proba(0.5,n) for n in Ln]
##plt.plot(Ln,Y)
##plt.show()