Programme modélisant salinisation d'une nappe phréatique

Version imprimable

Voici mon programme qui renvoie le temps à partir duquel la nappe est salinisée. Je voudrais aussi remplir graphiquement mon aquifère au fur et à mesure que le programme tourne, jusqu'au puit, ce qui modéliserait l'avancée de l'eau salée jusqu'au puits de pompage. J'ai essayé de faire cela en mettant :

Code:

1
2
3
4
5
   L=[ii for ii in range(n,Eau_salee+1,1)]   
 
    a1=[jj for jj in range(Eau_salee,n,1)]+L[::-1]
    b1=Si+Sj
    plt.fill_between(o,a1,b1,where=None,interpolate=True,hold=None)

Cependant voici l'erreur que python me renvoie : Traceback (most recent call last):

Code:

1
2
3
4
5
6
7
8
9
 File "<stdin>", line 1, in <module>
  File "<stdin>", line 46, in temps_salinisation
  File "C:\Users\Roxane\Anaconda\lib\site-packages\matplotlib\pyplot.py", line 2820, in fill_between
    interpolate=interpolate, **kwargs)
  File "C:\Users\Roxane\Anaconda\lib\site-packages\matplotlib\axes\_axes.py", line 4318, in fill_between
    y2 = ma.masked_invalid(self.convert_yunits(y2))
  File "C:\Users\Roxane\Anaconda\lib\site-packages\numpy\ma\core.py", line 2237, in masked_invalid
    condition = ~(np.isfinite(a))
TypeError: ufunc 'isfinite' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

En revanche, cela fonctionne en scilab en mettant :

Code:

1
2
3
        a1=[Leau_salee:n, n:-1:Leau_salee+1]
        b1=[Si,Sj]
        xfpoly(a1,b1,22)

mais, moi je dois le traduire en python.

Voici mon programme complet :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.pyplot import fill_between
 
x=[58,52,51,49,47,45,44,42,35,26,43,17,0]
y=[-50,-62,-61,-60,-60,-61,-61,-61,-61,-61,-61,-62,-61]
z=[-80,-88,-85,-83,-82,-83,-85,-85,-86,-87,-89,-91,-93]
a=[-89,-93,-95,-101,-112,-120,-126,-127,-124,-119,-118,-118,-117]
b=[-145,-155,-160,-168,-182,-192,-201,-204,-206,-207,-210,-212,-214]
c=[-379,-390,-400,-395,-396,-382,-373,-383,-434,-490,-537,-409,-409]
e=[-571,-549,-530,-519,-509,-495,-482,-494,-549,-611,-670,-504,-504]
 
o=[0,1,2,3,4,5,6,7,8,9,10,11,12]
 
px = [6,6]   
py= [0,-100] 
 
def temps_salinisation(n,q,p): #n distance du forage à la mer, q point du forage 
    plt.plot(x)
    plt.plot(y)
    plt.plot(z)
    plt.plot(a)
    plt.plot(b)
    plt.plot(c)
    plt.plot(e)
    px = [10,10] #position du puit
    py= [0,-120] #position du puit
    plt.plot(px,py, 'r-')
    plt.grid(True)
    plt.show()
 
    K=0.0006 #coefficient de perméabilité du sable fin 
    I=Y[p]
    J=Z[p]
    Eau_salee=n
    les=n #variable qui détermine les dx parcouru par l'eau à chaque dt jusqu'à atteindre Leau_salle-1. à t=0, les=Leau_salee
    Ymin=I[0]
    Zmax=J[0]
    for k in range(2,n,1):
        if I[k]< Ymin: 
            Ymin=I[k] #on détermine la profondeur maximale du toit de la nappe
        elif J[k]>Zmax:
            Zmax=J[k] # on détermine la hauteur maximale du plancher de la nappe
    hauteur_forage=[X[p][q]+.05,(Ymin+Zmax)/2] # détermine la profondeur du forage
    h=abs(Ymin-hauteur_forage[1]) #détermine la hauteur du cylindre dans lequel s'applique la loi de Darcy, en prenant en compte Ymin et Xmaxdétermine la profondeur du forage
    d= n-q #distance à parcourir jusqu'au point de forage.
    t=0
    while Eau_salee > q: #tant que la limite d'eau salée n'a pas atteint le puit
        Si=[]
        Sj=[]
        while les > Eau_salee-1 : 
            les= les-K*(h/d)*k #La loi de Darcy nous donne un débit volumique. on peut donc en déduire le volume d'eau introduit à chaque dt, et en simplifiant par la Section, on trouve la distance parcourue -tps nous donne le nombre de dt écoulés pour que l'eau salée atteigne le forage
            t=t+1
        Eau_salee=Eau_salee-1
        for i in range(Eau_salee,n,1) : 
            Sj=Sj+[J[n-i+1]]
        for i in range(Eau_salee,n,1): 
             Si=Si+[I[i]]   
    L=[ii for ii in range(n,Eau_salee+1,1)]   
 
    a1=[jj for jj in range(Eau_salee,n,1)]+L[::-1]
    b1=Si+Sj
    plt.fill_between(o,a1,b1,where=None,interpolate=True,hold=None)
    return (t)

Je m'adresse ici pour savoir si quelqu'un trouverait l'origine de l'erreur, et s'il aurait une idée pour la corriger, merci d'avance. (Je suis une élève de prepa bcpst, et je fais ce programme dans le cadre du TIPE)

26/03/2015, 12h17
BufferBob

:tagcode:

Salut,

Ton code est illisible sans balises.

On y trouve tout de même une chose étrange:

Code:

1
2
3
4
 
Sj = []
for i in range(Eau_salee,n,1) :
    Sj = [Sj, J[n-i+1]]

Si je le refais:

Code:

1
2
3
4
5
6
7
8
 
>>> Sj = []
>>> J = ['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']
>>> for i in range(10):
...     Sj = [Sj, J[i]]
... 
>>> Sj
[[[[[[[[[[[], 'f'], 'h'], 'j'], 'k'], 'l'], 'p'], 'p'], 'i'], 'g'], 'f']

c'est vraiment ça que tu veux ?

26/03/2015, 20h22
timshel

Bonsoir, alors tout d'abord, merci d'avoir répondu aussi vite.
Ensuite, tu as effectivement mis le doigt sur le problème, bien joué. Alors non effectivement, ce n'est pas ce que je veux. En suivant ton exemple, je voudrais qu'il me renvoie : ['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']. Comment modifierais tu ton exemple pour qu'il renvoie cela?
J'avais pensé à

Code:

Sj = Sj+[J[i]]]

mais cela ne fonctionne pas non plus...
27/03/2015, 04h44
papyfouette
Salut à toi,

Pourquoi tu ne tentes pas un simple Sj.append(i) ?
Où alors j'ai rien compris ! :P
Code:

1 2 3 4 for i in J(Vu que range(10) = J ?) Sj.append(i) >>> Sj ['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']

Salut,

A priori, pas besoin de faire une boucle...

Code:

1
2
3
4
5
6
7
8
>>> J = ['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']
>>> Sj = J[:]
>>> Sj
['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']
>>> Sj = J
>>> Sj
['f', 'h', 'j', 'k', 'l', 'p', 'p', 'i', 'g', 'f']
>>>

- W

Salut, merci pour vos réponses, j'ai finalement mis cela :

Code:

1
2
3
4
for i in range(Eau_salee,n,1) : 
            Sj=[J[n-i+1]]
          for i in range(Eau_salee,n,1): 
            Si=[I[i]]

Si je ne mets pas le plt.fill_between, le programme fonctionne, il me renvoie bien un temps, en revanche lorsque je rajoute :

Code:

1
2
3
4
5
L=[ii for ii in range(n,Eau_salee+1,1)]   
 
          a1=[jj for jj in range(Eau_salee,n,1)]+L[::-1]
          b1=Si+Sj
          plt.fill_between(o,a1,b1,where=None,interpolate=True,hold=None)

qui est dans le but de remplir l'espace entre 2 courbes sur un graphe au fur et à mesure que le programme tourne, ce qui modélise en fait l'entrée du sel dans la nappe phréatique, et bien le programme ne marche plus et l'erreur est :

Code:

1
2
3
4
5
6
7
8
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "<stdin>", line 46, in temps_salinisation
  File "C:\Users\Roxane\Anaconda\lib\site-packages\matplotlib\pyplot.py", line 2820, in fill_between
    interpolate=interpolate, **kwargs)
  File "C:\Users\Roxane\Anaconda\lib\site-packages\matplotlib\axes\_axes.py", line 4331, in fill_between
    raise ValueError("Argument dimensions are incompatible")
ValueError: Argument dimensions are incompatible

Je ne sais pas comment faire pour avoir les bonnes dimensions, si ça trouve je ne peux pas faire ce que je veux avec fill_betwenn, mais alors auriez vous une autre idée pour remplir l'espace entre deux courbes au fur et à mesure que mon programme tourne ?

Je vous mets le programme en entier avec les modifications apportées :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.pyplot import fill_between
 
X=[[58,52,51,49,47,45,44,42,35,26,43,17,0],[58,50,40,37,37,35,36,35,28,26,38,8,0],[56,51,41,30,27,25,24,24,19,21,22,4,0],[49,46,40,29,20,15,12,12,10,18,11,3,0]]
Y=[[-50,-62,-61,-60,-60,-61,-61,-61,-61,-61,-61,-62,-61],[-47,-60,-63,-62,-60,-60,-61,-61,-60,-60,-61,-60,-59],[-45,-55,-62,-63,-54,-60,-60,-60,-60,-61,-59,-57,-55],[-43,-49,-57,-61,-62,-61,-61,-60,-60,-58,-56,-53,-52]]
Z=[[-80,-88,-85,-83,-82,-83,-85,-85,-86,-87,-89,-91,-93],[-72,-81,-83,-82,-79,-70,-82,-83,-83,-85,-87,-88,-88],[-65,-72,-78,-79,-70,-77,-78,-80,-82,-84,-84,-84,-84],[-57,-63,-70,-74,-75,-74,-76,-77,-79,-80,-79,-78,-79]]
A=[[-89,-93,-95,-101,-112,-120,-126,-127,-124,-119,-118,-118,-117],[-84,-88,-91,-95,-103,-104,-124,-125,-122,-117,-116,-115,-115],[-74,-83,-92,-92,-87,-108,-118,-123,-120,-116,-114,-113,-112],[-62,-85,-93,-92,-92,-97,-110,-118,-118,-116,-112,-109,-108]]
B=[[-145,-155,-160,-168,-182,-192,-201,-204,-206,-207,-210,-212,-214],[-135,-147,-155,-162,-173,-177,-199,-204,-206,-206,-209,-210,-212],[-123,-140,-154,-159,-157,-181,-194,-203,-205,-206,-208,-209,-209],[-110,-140,-155,-158,-162,-171,-188,-199,-205,-207,-206,-206,-207]]
C=[[-379,-390,-400,-395,-396,-382,-373,-383,-434,-490,-537,-409,-409],[-335,-342,-348,-352,-343,-338,-347,-369,-417,-471,-416,-416,-416],[-281,-299,-313,-316,-313,-329,-335,-359,-405,-455,-419,-419,-419],[-209,-214,-230,-256,-282,-287,-292,-320,-370,-420,-420,-420,-420]]
D=[[-571,-549,-530,-519,-509,-495,-482,-494,-549,-611,-670,-504,-504],[-541,-523,-493,-478,-468,-444,-461,-482,-532,-592,-515,-515,-515],[-507,-483,-470,-471,-441,-449,-449,-474,-519,-573,-519,-519,-519],[-589,-504,-459,-453,-447,-426,-405,-420,-470,-520,-520,-520,-520]]
 
x=[58,52,51,49,47,45,44,42,35,26,43,17,0]
y=[-50,-62,-61,-60,-60,-61,-61,-61,-61,-61,-61,-62,-61]
z=[-80,-88,-85,-83,-82,-83,-85,-85,-86,-87,-89,-91,-93]
a=[-89,-93,-95,-101,-112,-120,-126,-127,-124,-119,-118,-118,-117]
b=[-145,-155,-160,-168,-182,-192,-201,-204,-206,-207,-210,-212,-214]
c=[-379,-390,-400,-395,-396,-382,-373,-383,-434,-490,-537,-409,-409]
e=[-571,-549,-530,-519,-509,-495,-482,-494,-549,-611,-670,-504,-504]
 
o=[0,1,2,3,4,5,6,7,8,9,10,11,12]
 
px = [6,6]   
py= [0,-100] 
 
def temps_salinisation(n,q,p): #n distance du forage à la mer, q point du forage 
    plt.plot(x)
    plt.plot(y)
    plt.plot(z)
    plt.plot(a)
    plt.plot(b)
    plt.plot(c)
    plt.plot(e)
    px = [10,10] #position du puit
    py= [0,-120] #position du puit
    plt.plot(px,py, 'r-')
    plt.grid(True)
    plt.show()
 
    K=0.0006 #coefficient de perméabilité du sable fin 
    I=Y[p]
    J=Z[p]
    Eau_salee=n
    les=n #ensemble des dx jusqu'à ce qu'on atteigne 1 
    Ymin=I[0]
    Zmax=J[0]
    for k in range(2,n,1):
        if I[k]< Ymin: 
            Ymin=I[k] #on détermine la profondeur maximale du toit de la nappe
        elif J[k]>Zmax:
            Zmax=J[k] # on détermine la hauteur maximale du plancher de la nappe
    hauteur_forage=[X[p][q]+.05,(Ymin+Zmax)/2] # détermine la profondeur du forage
    h=abs(Ymin-hauteur_forage[1]) #détermine la hauteur du cylindre dans lequel s'applique la loi de Darcy, en prenant en compte Ymin et Xmaxdétermine la profondeur du forage
    d= n-q #distance à parcourir jusqu'au point de forage.
    t=0
    while Eau_salee > q: #tant que la limite d'eau salée n'a pas atteint le puit
        Si=[]
        Sj=[]
        while les > Eau_salee-1 : 
            les= les-K*(h/d)*k #La loi de Darcy nous donne un débit volumique. on peut donc en déduire le volume d'eau introduit à chaque dt, et en simplifiant par la Section, on trouve la distance parcourue -tps nous donne le nombre de dt écoulés pour que l'eau salée atteigne le forage
            t=t+1
        Eau_salee=Eau_salee-1
        for i in range(Eau_salee,n,1) : 
            Sj=[J[n-i+1]]
        for k in range(Eau_salee,n,1): 
            Si=[I[k]]  
    L=[ii for ii in range(n,Eau_salee+1,1)]   
 
    a1=[jj for jj in range(Eau_salee,n,1)]+L[::-1]
    b1=Si+Sj
    plt.fill_between(o,a1,b1,where=None,interpolate=True,hold=None)
    return (t)

28/03/2015, 12h33
marco056

Ton programme ne retourne pas d'erreur en python3.
28/03/2015, 14h03
timshel

Ah oui ?? Et on a un remplissage qui se fait au fur et à mesure que le programme tourne? Ou seulement les courbes s'affichent et aucun remplissage n'est visible?
29/03/2015, 19h18
marco056

Citation:

Envoyé par timshel

Ah oui ?? Et on a un remplissage qui se fait au fur et à mesure que le programme tourne? Ou seulement les courbes s'affichent et aucun remplissage n'est visible?

Aucun affichage, mais pas d'erreur de déclenchée.