sys-bio / tellurium

Python Environment for Modeling and Simulating Biological Systems
http://tellurium.analogmachine.org/
Apache License 2.0
106 stars 36 forks source link

Trying to pretty print r.values() results in runtime error #357

Open eugenio opened 6 years ago

eugenio commented 6 years ago

Not sure if this is the right place to post this issue, I case please remove it and i will reopen it in roadrunner...

this line

pp.pprint(r.values())

of this file

import matplotlib.pyplot as plt
import pprint as pp
# import numpy as np
import tellurium as te

te.setDefaultPlottingEngine('matplotlib')

axarr = plt.subplots(3, 1)
r = te.loada("model_ComplexI.txt")

# r.model['init(Q_Rd)'] = 1 * 10 ^ -6
# r.model['init(Q_Ox)'] = 1 * 10 ^ -6
# r.model['init(NAD_Rd)'] = 1 * 10 ^ -3
# r.model['init(NAD_Ox)'] = 1 * 10 ^ -3
print("I_nc = " + str(r.getValue("I_nc")))
print("r_kiqox_Qox = " + str(r.getValue("r_kiqox_Qox")))
print("r_kiqrd_Qrd = " + str(r.getValue("r_kiqrd_Qrd")))
print(r.J0, r.J1)
print("E_cell = " + str(r.getValue("E_cell")))
print("E_pmf = " + str(r.getValue("E_pmf")))
print("Q_reaction = " + str(r.getValue("Q_reaction")))

pp.pprint(r.values())
# r.getIntegratorByName('gillespie').seed = 12345
m = r.simulate(0, 100, 1000)
# n = r.simulate( 0, 0, 1000, selections=["Time", "NAD_Rd"])
# o = r.simulate( 0, 0, 1000, selections=["Time", "Q_Ox"])

plt.subplot(3,1,1)
m_plot = plt.plot(m[:,1])
# plt.subplot (3,1,2)
# n_plot = plt.plot(n[:,1])
# plt.subplot (3,1,3)
# o_plot = plt.plot(o[:,1])

plt.show()

results in :

Traceback (most recent call last):
  File "/home/*********/***********/synthetic biology/***************************/ComplexI-sim.py", line 24, in <module>
    pp.pprint(r.values())
  File "/home/*********/***********/miniconda3/envs/tellurium/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 3533, in values
    return [self.getValue(k) for k in self.keys(types)]
  File "/home/*********/***********/miniconda3/envs/tellurium/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 3533, in <listcomp>
    return [self.getValue(k) for k in self.keys(types)]
  File "/home/*********/***********/miniconda3/envs/tellurium/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 3379, in getValue
    return _roadrunner.RoadRunner_getValue(self, *args)
RuntimeError: Invalid id 'n' for floating initial amount```

Model is: 

model complexI();

    compartment cytoplasm = 1 , membrane = 0.1, extracellular = 10^6 litres;

    species Enz, NAD_Ox, NAD_Rd, Q_Ox, Q_Rd, Hrea;

    //reaction 
    J0: NAD_Ox + Enz -> NAD_Rd + Enz ; -Kvel*Enz*NAD_Ox;
    J1: Q_Rd + Enz -> Q_Ox + Enz; -Kvel*Enz*Q_Rd ;
    J2: Hext -> 4.5Hrea; -Kvel*Enz;

    Enz in cytoplasm;
    Hcyt in cytoplasm;
    Hext in extracellular;
    Q_Rd in membrane;
    NAD_Rd in cytoplasm;
    Q_Ox in membrane;
    NAD_Ox in cytoplasm;
    const   Hrea in cytoplasm;

    Q_Rd = 1 * 10 ^ -6 mole;
    Q_Ox = 1 * 10 ^ -6 mole;
    NAD_Rd = 1 * 10 ^ -3 mole;
    NAD_Ox = 1 * 10 ^ -3 mole;

    n = 2;
    F = 96485.3329;
    R = 8.314472;
    T = 310.15;

    Enz = 0.0000000017 mole;
    Hcyt := 10^-7.5 +(-Hrea); 
    Hext = 10^-3.5 M;
    Hrea = 0.0 mole

    unit Krate = mole/second;
    k_Ox = 0.014859782 Krate;
    k_Rd = 0.000014859782 Krate;
    k_NAD_Rd = 6.1 microM;
    k_Q_Ox = 13 microM;
    ki_Q_Ox = 1367 microM;
    ki_Q_Rd = 56.8 microM;
    k_NAD_Ox = 2064 microM;
    k_Q_Rd = 4 microM;

    E0qred_qox = 0.045 V ;
    E0nadr_nadox = -0.320 V ; 

    E0cell = -(E0qred_qox) + (E0nadr_nadox);

    k_Rd_neq := 1.54 *10^+11 / ((k_NAD_Rd / k_NAD_Ox) * (k_Q_Rd/k_Q_Rd) * k_Ox);
    k_Ox_neq := 1.54 *10^+11 / ((k_NAD_Rd / k_NAD_Ox) * (k_Q_Rd/k_Q_Rd) * k_Rd);
    Delta_pH_o_i := -log10(Hext)-(-log10(Hcyt));
    E_pmf := (2.303*R*T*Delta_pH_o_i)/F ;
    Q_reaction := (NAD_Ox / NAD_Rd) * ( Q_Rd / Q_Ox) ;
    E_cell := E0cell - ((((2.303*R*T)/F)/n) * log10(Q_reaction)) - E_pmf ;
    pH_cyt := -log10(Hcyt) ;

    r_kiqox_Qox := Q_Ox/ki_Q_Ox
    r_kiqrd_Qrd := Q_Rd/ki_Q_Rd

    I_nc := (1/(1 +( r_kiqox_Qox) + ( r_kiqrd_Qrd)))^2 ;
    Kvel := E_cell* I_nc *(((k_Ox_neq*(NAD_Rd/k_NAD_Rd)*(Q_Ox/k_Q_Ox))-(k_Rd_neq * (NAD_Ox / k_NAD_Ox )*( Q_Rd / k_Q_Rd )))/((1+(NAD_Rd/k_NAD_Rd)+ (NAD_Ox / k_NAD_Ox))*(1+(Q_Ox/ k_Q_Ox)+(Q_Rd/k_Q_Rd))));

end

If I remove the said line or comment it out the above script result in : 

RuntimeError: CVODE Error: CV_CONV_FAILURE: Convergence test failures occurred too many times (= MXNCF = 10) during one internal timestep or occurred with |h| = hmin.; In virtual double rr::CVODEIntegrator::integrate(double, double)



Probably my model has other issues but I don't understand the error message and what can I do to fix my model eventually
clayton-r commented 5 years ago

I have exactly the same issue, did you ever solve it?

hsauro commented 5 years ago

What is the actual problem? If I type r.values() on a model I get a listing of values, after which I can still run a simulation. Do you have a small example that shows the issue?

On Mon, Aug 19, 2019 at 4:19 PM Zer0Credibility notifications@github.com wrote:

I have exactly the same issue, did you ever solve it?

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/sys-bio/tellurium/issues/357?email_source=notifications&email_token=AAIBSDSA6J4Y7T2LN2MJLV3QFMS7BA5CNFSM4EYGGRI2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD4USJ4A#issuecomment-522790128, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIBSDS2F4R5LXKQNGXZRG3QFMS7BANCNFSM4EYGGRIQ .

-- Herbert Sauro, Associate Professor University of Washington, Bioengineering 206-685-2119, www.sys-bio.org hsauro@uw.edu Books: http://books.analogmachine.org/

clayton-r commented 5 years ago

Here's an example of a run that reproduces this behavior:

model: https://raw.githubusercontent.com/Zer0Credibility/Rhodo/master/CGA009.an

r = te.tellurium.loadAntimonyModel('CGA009.an')
m = r.simulate(0, 300, 300)

And the output:

Error: CVODE Error: CV_CONV_FAILURE, Module: CVODE, Function: CVode, Message: At t = 0 and h = 2.20348e-46, the corrector convergence test failed repeatedly or with |h| = hmin.
Error: rr::IntegratorException::IntegratorException(const std::string &, const std::string &)what: CVODE Error: CV_CONV_FAILURE: Convergence test failures occurred too many times (= MXNCF = 10) during one internal timestep or occurred with |h| = hmin., where: virtual double rr::CVODEIntegrator::integrate(double, double)
Traceback (most recent call last):
  File "/Users/User/Documents/MSA/Simulate.py", line 6, in <module>

  File "/Users/User/anaconda3/envs/SBML36/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 2913, in simulate
    result = self._simulate(o)
  File "/Users/User/anaconda3/envs/SBML36/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 2628, in _simulate
    def _simulate(self, *args): return _roadrunner.RoadRunner__simulate(self, *args)
RuntimeError: CVODE Error: CV_CONV_FAILURE: Convergence test failures occurred too many times (= MXNCF = 10) during one internal timestep or occurred with |h| = hmin.; In virtual double rr::CVODEIntegrator::integrate(double, double)

Process finished with exit code 1
hsauro commented 5 years ago

We'll check it out. It's a big one. Do you know if anyone has managed to simulate this on any platform?

Herbert

On Mon, Aug 26, 2019 at 4:21 PM Zer0Credibility notifications@github.com wrote:

Here's an example of a run that reproduces this behavior:

model: https://raw.githubusercontent.com/Zer0Credibility/Rhodo/master/CGA009.an

r = te.tellurium.loadAntimonyModel('CGA009.an') m = r.simulate(0, 300, 300)

And the output:

Error: CVODE Error: CV_CONV_FAILURE, Module: CVODE, Function: CVode, Message: At t = 0 and h = 2.20348e-46, the corrector convergence test failed repeatedly or with |h| = hmin. Error: rr::IntegratorException::IntegratorException(const std::string &, const std::string &)what: CVODE Error: CV_CONV_FAILURE: Convergence test failures occurred too many times (= MXNCF = 10) during one internal timestep or occurred with |h| = hmin., where: virtual double rr::CVODEIntegrator::integrate(double, double) Traceback (most recent call last): File "/Users/User/Documents/MSA/Simulate.py", line 6, in

File "/Users/User/anaconda3/envs/SBML36/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 2913, in simulate result = self._simulate(o) File "/Users/User/anaconda3/envs/SBML36/lib/python3.6/site-packages/roadrunner/roadrunner.py", line 2628, in _simulate def _simulate(self, args): return _roadrunner.RoadRunner__simulate(self, args) RuntimeError: CVODE Error: CV_CONV_FAILURE: Convergence test failures occurred too many times (= MXNCF = 10) during one internal timestep or occurred with |h| = hmin.; In virtual double rr::CVODEIntegrator::integrate(double, double)

Process finished with exit code 1

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/sys-bio/tellurium/issues/357?email_source=notifications&email_token=AAIBSDSCGWSQB44D6A2NJRTQGRQQNA5CNFSM4EYGGRI2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD5F7B5Q#issuecomment-525070582, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIBSDTQT2PRHSW5HNXDXW3QGRQQNANCNFSM4EYGGRIQ .

-- Herbert Sauro, Associate Professor University of Washington, Bioengineering 206-685-2119, www.sys-bio.org hsauro@uw.edu Books: http://books.analogmachine.org/

clayton-r commented 5 years ago

Thanks for having a look! I've managed to get this model running on Cobra -- I've also found that I get the same error on both Ubuntu 18.04 and MacOS Mojave

hsauro commented 5 years ago

Does cobra do kinetic simulation?

Herbert

On Mon, Aug 26, 2019 at 6:04 PM Zer0Credibility notifications@github.com wrote:

Thanks for having a look! I've managed to get this model running on Cobra -- I've also found that I get the same error on both Ubuntu 18.04 and MacOS Mojave

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/sys-bio/tellurium/issues/357?email_source=notifications&email_token=AAIBSDQ2UGRJ6K4TD56COWDQGR4PPA5CNFSM4EYGGRI2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD5GEG3A#issuecomment-525091692, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIBSDVWXPHOPSOIMVBP2JTQGR4PPANCNFSM4EYGGRIQ .

-- Herbert Sauro, Associate Professor University of Washington, Bioengineering 206-685-2119, www.sys-bio.org hsauro@uw.edu Books: http://books.analogmachine.org/

clayton-r commented 5 years ago

No I don't think so, at least not currently.

It looks like setting the number of allowed convergence failures might help.

CVodeSetMaxConvFails
hsauro commented 5 years ago

Ok I think I understand, you manage to carry out a flux balance analysis on cobra?

Where did the kinetic rate laws come from in the model?

Herbert

On Mon, Aug 26, 2019 at 6:46 PM Zer0Credibility notifications@github.com wrote:

No I don't think so, at least not currently.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/sys-bio/tellurium/issues/357?email_source=notifications&email_token=AAIBSDROAQDF7MIQMPYG7ALQGSBPRA5CNFSM4EYGGRI2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD5GGJ7A#issuecomment-525100284, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIBSDXYT7RXF3WQZOPXPELQGSBPRANCNFSM4EYGGRIQ .

-- Herbert Sauro, Associate Professor University of Washington, Bioengineering 206-685-2119, www.sys-bio.org hsauro@uw.edu Books: http://books.analogmachine.org/

clayton-r commented 5 years ago

Yeah Flux balance works with the model with some slight adjustments to add objective coefficients, etc.

I used parameter balancing to generate the kinetic components.

hsauro commented 5 years ago

Ok I see.

There could be all sorts of problems that might result in the problems you found. What I would do is break the model into more manageable sub models (less than 50 reactions on each) and make sure each sub part can sustain a steady state. We’ll have a look at our end but I suspect you have a model that doesn’t have a solution. The first thing to look at are the reaction rates and rates of change at time zero. You can check whether the rates of change or reaction rates with values of infinity. If these are present the integrator will fail.

Herbert

On Mon, Aug 26, 2019 at 6:59 PM Zer0Credibility notifications@github.com wrote:

Yeah Flux balance works with the model with some slight adjustments to add objective coefficients, etc.

I used parameter balancing to generate the kinetic components.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/sys-bio/tellurium/issues/357?email_source=notifications&email_token=AAIBSDSJOXH3ILM66K33NR3QGSC6DA5CNFSM4EYGGRI2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD5GG57A#issuecomment-525102844, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIBSDTLHZ4UHMRQW444KRDQGSC6DANCNFSM4EYGGRIQ .

-- Herbert Sauro, Associate Professor University of Washington, Bioengineering 206-685-2119, www.sys-bio.org hsauro@uw.edu Books: http://books.analogmachine.org/