Integration von functools.partial Objekte

Question

Ich versuche, einen stochastischen Simulator zu machen. Allerdings möchte ich in einer deterministischen Umgebung einige Plausibilitätsprüfungen durchführen, daher muss ein System von ODEs integriert werden. Speziell das Produkt von self.stoich_mat und prop_cont wie in der integrand Funktion angegeben. Aus irgendeinem Grund die Integra Funktion arbeitet, Ausgabe geben [15. -15.], aber scipy's integrate.odeint nicht und den Fehler wirft mich nicht funktioniert:

[ 15. -15.] #THIS IS THE OUTPUT FROM THE INTEGRAND FUNCTION 

Traceback (most recent call last): 
    File "C:/Users/Tomek/Documents/MasterThesis/sim_method5.py", line 474, in <module> 
    Simulator2.run() 
    File "C:/Users/Tomek/Documents/MasterThesis/sim_method5.py", line 455, in run 
    full_output=True) 
    File "C:\ProgramData\Anaconda2\lib\site-packages\scipy\integrate\odepack.py", line 215, in odeint 
    ixpr, mxstep, mxhnil, mxordn, mxords) 
TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'safe' 

ich die MarkoV Klasse gebaut mit verschiedenen Simulationstypen zu arbeiten, und deshalb Die self.propensity_generator könnte ein wenig kompliziert aussehen. Der folgende Code ist der kleinste Betrag, der immer noch den Fehler anzeigt. Ich habe versucht, Lambda-Ausdrücke anstelle von partiellen functools zu verwenden, aber das scheint nicht zu helfen.

Der Code:

class MarkoV(object): 
    def __init__(self,rates,stoich_mat): 
     self.stoich_mat=stoich_mat 
     self.rates=rates 
     self.nr_reactions=rates.shape[1] 

    def init(self,init_state,tf,t0=0.0): 
     self.state=init_state 
     self.init_state=np.copy(init_state) 
     self.nr_species=init_state.shape[0] 
     self.tf=tf 
     self.t_init=t0 
     self.propensity_generator() 

    def uni_het_bi_mol_const(self,state,t,k,el): 
     return k*np.prod(state[el]) 
    def uni_het_bi_mol_cat_t(self,state,t,k,cat,el): 
     return k*np.prod(state[el])*self.init_state[cat][0](t=t) 

    def propensity_generator(self): 
     self.propensity_funcs=range(self.stoich_mat.shape[1]) 
     for i,reac in enumerate(self.stoich_mat.T): 
      j=np.where(reac<0)[0] 
      if isinstance(self.rates[0, i], float): 
       k=self.rates[0, i] 
       self.propensity_funcs[i] = functools.partial(self.uni_het_bi_mol_const,k=k,el=j) 
      else: 
       k,cat=self.rates[0,i] 
       #HERE THE DETERMINISTIC INPUT IS DEFINED, ESSENTIAL TO MY PROBLEM, BUT MAYBE A HURDLE 
       name_func,func_kwargs=self.init_state[cat][0] 
       self.init_state[cat][0]=functools.partial(getattr(self,name_func),**func_kwargs) 
       #self.propensity_funcs[i]=lambda state,t:self.uni_het_bi_mol_cat_t(state=state,t=t,k=k,cat=cat,el=j) #TRY LAMBDA 
       self.propensity_funcs[i] = functools.partial(self.uni_het_bi_mol_cat_t,k=k,cat=cat,el=j) 
     self.propensity_funcs=np.array(self.propensity_funcs) 

class deterministic(MarkoV): 
    def init_sym(self,init_state,tf,species_continuous,t0=0.0): 
     self.species_continuous=species_continuous 
     self.init(init_state,tf,t0) 
     self.species_discrete=[i for i in np.arange(self.nr_species) if i not in self.species_continuous] 

    def sinput(self,t,amplitude=6.0,frequency=0.05,offset=1.0): 
     return amplitude*np.sin(frequency*t)+amplitude+offset 

    def integrand(self,state,t): 
     prop_cont=np.zeros(self.nr_reactions) 
     for i,func in enumerate(self.propensity_funcs): 
      prop_cont[i]=func(state=state,t=t) 
     return np.dot(self.stoich_mat[self.species_discrete,:],prop_cont) 

    def run(self,t_step=1e-2): 
     t=np.arange(self.t_init,self.tf,t_step) 

     h=integrate.odeint(func=self.integrand, 
          y0=self.state[self.species_discrete], 
          t=t, 
          full_output=True) 

if __name__ == '__main__': 
    import scipy.integrate as integrate 
    import numpy as np 
    import functools 
    import numbers 

    rates=np.array([(1.2,0),6.0],dtype=object,ndmin=2) 
    stoich_mat=np.array([[ 0, 0], # input 
         [-1, 1], # A 
         [ 1,-1]])# A* 
    init_state=np.array([("sinput",{"amplitude":2.0,"frequency":0.05,"offset":0.5}),5.0,5.0],ndmin=2).T 
    species_continuous=[0] 

    Simulator2=deterministic(rates=rates,stoich_mat=stoich_mat) 
    Simulator2.init_sym(init_state=init_state,tf=20.0,species_continuous=species_continuous) 
    print Simulator2.integrand(state=init_state,t=0) 
    print "" 
    Simulator2.run() 

Answer 1

In [126]: init_state=np.array([("sinput",{"amplitude":2.0,"frequency":0.05,"offs 
    ...: et":0.5}),5.0,5.0],ndmin=2).T 
In [127]: init_state 
Out[127]: 
array([[('sinput', {'frequency': 0.05, 'amplitude': 2.0, 'offset': 0.5})], 
     [5.0], 
     [5.0]], dtype=object) 
In [128]: init_state.shape 
Out[128]: (3, 1) 
In [129]: species_discrete=[i for i in np.arange(3) if i not in [0]] 
In [130]: species_discrete 
Out[130]: [1, 2] 
In [131]: init_state[species_discrete] 
Out[131]: 
array([[5.0], 
     [5.0]], dtype=object) 
In [132]: _.shape 
Out[132]: (2, 1) 

so gut ich kann Ihnen sagen, setzen y0 in odeint auf eine (2,1) Objekt-Array.

ich vermuten, dass Ihr Code laufen würde, oder zumindest nicht diesen Fehler erhöhen, wenn es waren:

In [133]: init_state[species_discrete].astype(float).ravel() 
Out[133]: array([ 5., 5.]) 

Sie auch sicher, dass

Simulator2.integrand(y0, 0) 

Läufe machen wollen, und gibt ein ähnliches (2,) Array.