scipy.integrate.ode et numpy.view
Bonjour à tous
Je cherche à résoudre 2 systèmes d'équations couplées. Appelons les système A et système B.
L'un de ces 2 systèmes (disons système A) est un système EDO. Pour éviter de dupliquer l'information commune au 2 systèmes, je souhaite procéder avec des numpy.view
Je m'explique via le code :
Code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
| import numpy as np
import scipy
t0,t1,dt = 0.0,5.0, 1.0
data = np.ones((5,2))
data[:,1]*=2
y=np.array([0.0,0.0]) ### no matter default value
r = scipy.integrate.ode(f)
r.set_integrator('dopri5', rtol=1e-3, atol=1e-6 )
r.set_f_params(0.05)
#r.set_initial_value(y, t0); r._y = data[2] ### Apparently equivalent
r.set_initial_value(data[2], t0) ### Apparently equivalent
print(np.shares_memory(r.y,y))
print(np.shares_memory(r.y,data)) |
Ici, à l'état initial, j'ai donc une synchronisation entre r.y (le système A) et data[2] (data contenant les informations du système B). Si je modifie l'un, l'autre se modifie aussi et vis versa. D'ailleurs, la commande r.y.base nous montre bien que r.y n'est qu'une view de data.
Ce comportement m'arrange bien car si le système B auquel le système A est couplé, vient modifier une valeur, elle est directement mis à jour.
Jusque là tout va bien. Maintenant, je fais avancer mon EDO dans le temps
Code:
1
2
3
4
5
6
| while r.successful() and r.t < t1:
r.integrate(r.t+dt, step=True)
print(r.t+dt,r.y)
print(np.shares_memory(r.y,data))
print(data) |
et là drame, data et r.y ne sont plus synchro. r.y n'est plus une view de data.
Il semblerait que la fonction integrate créée une nouvelle instance de son attribut r.y plutot que de simplement le mettre à jour. J'ai été voir un peu le code de cette fonction
Code:
https://github.com/scipy/scipy/blob/v0.19.1/scipy/integrate/_ode.py#L396
mais en remontant un peu ça renvoie a du code fortran et là j'y comprends plus rien.
Comment puis-je résoudre (ou contourner) ce problème autrement qu'en me tapant une copie de l'information entre r.y et data, dont je devrais gérer la synchronisation à la main ?
Merci
