Pysb Import (#463)

* feature(core) allow observables to have expressions w * feature(python) add native pysb import and add tests * fix(tests) update pregenerated models * fix(tests) update tests * fix(python,sbml_import) fix wrong jacobian, it always suprises me what you can get away with ... * fix(tests) bugfixes after tests and install script for bngl * fix(tests) update robertson example * fix(tests) update neuron models * fix(doc) update documentation * fix(style) fix codacy issues * fix(style) fix undefined import error * fix(python) move code to pysb_imort, split into multiple functions * remove N_Vector variants of fw, fdwdp, fdwdx * fix(ode_export) fix problems for nw==0 and improve assignment of numeric values to ModelQuantity values * feature(tests) add more testcases for pysb import * fix(tests) readd compare_to_sbml_import testcase * fix(style) fix whitespace formatting in testPYSB * fix(python) fix unused import * increase version number
AMICI-dev · Nov 17, 2018 · 11c8ec3 · 11c8ec3
1 parent ef26147
commit 11c8ec3
Show file tree

Hide file tree

Showing 55 changed files with 863 additions and 330 deletions.
diff --git a/.gitignore b/.gitignore
@@ -193,8 +193,15 @@ python/test.txt
 tests/test/*
 tests/test_model_steadystate_scaled/*
 tests/test_model_presimulation/*
-
-
+tests/test_model_presimulation_pysb/*
+tests/test_model_presimulation_sbml/*
+tests/bax_pore_sequential/*
+tests/bax_pore/*
+tests/bngwiki_egfr_simple/*
+tests/earm_1_3/*
+tests/expression_observables/*
+tests/robertson/*
+tests/tyson_oscillator/*
 
 coverage_py.xml
 
@@ -217,3 +224,5 @@ testSuite.txt
 python/sdist/build/*
 *.so
 python/sdist/amici/git_version.txt
+
+
diff --git a/.travis.yml b/.travis.yml
@@ -55,6 +55,8 @@ install:
     - ./scripts/buildSuiteSparse.sh
     - ./scripts/buildSundials.sh
     - ./scripts/buildCpputest.sh
+    - ./scripts/buildBNGL.sh
+    - export BNGPATH=$(pwd)/ThirdParty/BioNetGen-2.3.2
     - ./scripts/buildAmici.sh
 
 script:  

diff --git a/AMICI_guide.pdf b/AMICI_guide.pdf
diff --git a/include/amici/model.h b/include/amici/model.h
@@ -771,27 +771,27 @@ namespace amici {
          * in natural order.
          */
         void requireSensitivitiesForAllParameters();
-
+        
         /**
          * @brief Recurring terms in xdot
          * @param t timepoint
-         * @param x Vector with the states
+         * @param x array with the states
          */
-        void fw(const realtype t, const N_Vector x);
-
+        void fw(const realtype t, const realtype *x);
+        
         /**
          * @brief Recurring terms in xdot, parameter derivative
          * @param t timepoint
-         * @param x Vector with the states
+         * @param x array with the states
          */
-        void fdwdp(const realtype t, const N_Vector x);
+        void fdwdp(const realtype t, const realtype *x);
 
         /**
          * @brief Recurring terms in xdot, state derivative
          * @param t timepoint
-         * @param x Vector with the states
+         * @param x array with the states
          */
-        void fdwdx(const realtype t, const N_Vector x);
+        void fdwdx(const realtype t, const realtype *x);
 
         /** residual function
          * @param it time index
@@ -1195,8 +1195,9 @@ namespace amici {
          * @param p parameter vector
          * @param k constant vector
          * @param h heavyside vector
+         * @param w repeating elements vector
          **/
-        virtual void fy(realtype *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+        virtual void fy(realtype *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w) {
             throw AmiException("Requested functionality is not supported as (%s) is not implemented for this model!",__func__);
         }
 
@@ -1208,8 +1209,10 @@ namespace amici {
          * @param k constant vector
          * @param h heavyside vector
          * @param ip parameter index w.r.t. which the derivative is requested
+         * @param w repeating elements vector
+         * @param dwdp Recurring terms in xdot, parameter derivative
          **/
-        virtual void fdydp(realtype *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) {
+        virtual void fdydp(realtype *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp) {
             throw AmiException("Requested functionality is not supported as (%s) is not implemented for this model!",__func__);
         }
 
@@ -1220,8 +1223,10 @@ namespace amici {
          * @param p parameter vector
          * @param k constant vector
          * @param h heavyside vector
+         * @param w repeating elements vector
+         * @param dwdx Recurring terms in xdot, state derivative
          **/
-        virtual void fdydx(realtype *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+        virtual void fdydx(realtype *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx) {
             throw AmiException("Requested functionality is not supported as (%s) is not implemented for this model!",__func__);
         }
 

diff --git a/matlab/@amifun/getArgs.m b/matlab/@amifun/getArgs.m
@@ -56,7 +56,7 @@
         case 'root'
             this.argstr = ['(realtype *root, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h' dx ')'];
         case 'y'
-            this.argstr = '(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)';
+            this.argstr = '(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w)';
         case 'z'
             this.argstr = '(double *z, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)';
         case 'rz'
@@ -66,9 +66,9 @@
         case 'srz'
             this.argstr = '(double *srz, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *sx, const int ip)';
         case 'dydp'
-            this.argstr = '(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip)';
+            this.argstr = '(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp)';
         case 'dydx'
-            this.argstr = '(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)';
+            this.argstr = '(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx)';
         case 'dzdp'
             this.argstr = '(double *dzdp, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip)';
         case 'dzdx'
@@ -127,4 +127,4 @@
             %nothing
     end
 
-end
+end
diff --git a/matlab/examples/amiExamples.m b/matlab/examples/amiExamples.m
@@ -11,4 +11,6 @@
 example_adjoint_hessian
 example_dirac_adjoint_hessVecProd
 example_jakstat_adjoint_hvp
+example_robertson
+example_neuron
 path(oldpath);
diff --git a/models/model_dirac/model_dirac.h b/models/model_dirac/model_dirac.h
@@ -1,6 +1,6 @@
 #ifndef _amici_model_dirac_h
 #define _amici_model_dirac_h
-/* Generated by amiwrap (R2017b) acb0375e4e0c93fb725d81b42e7e818c9ba7b65d */
+/* Generated by amiwrap (R2017b)  */
 #include <cmath>
 #include <memory>
 #include "amici/defines.h"
@@ -26,15 +26,15 @@ extern void dJydy_model_dirac(double *dJydy, const int iy, const realtype *p, co
 extern void deltasx_model_dirac(double *deltasx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const int ip, const int ie, const realtype *xdot, const realtype *xdot_old, const realtype *sx, const realtype *stau);
 extern void deltax_model_dirac(double *deltax, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ie, const realtype *xdot, const realtype *xdot_old);
 extern void dxdotdp_model_dirac(realtype *dxdotdp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp);
-extern void dydx_model_dirac(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
+extern void dydx_model_dirac(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx);
 extern void qBdot_model_dirac(realtype *qBdot, const int ip, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *xB, const realtype *w, const realtype *dwdp);
 extern void root_model_dirac(realtype *root, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
 extern void sigmay_model_dirac(double *sigmay, const realtype t, const realtype *p, const realtype *k);
 extern void stau_model_dirac(double *stau, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *sx, const int ip, const int ie);
 extern void sxdot_model_dirac(realtype *sxdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *sx, const realtype *w, const realtype *dwdx, const realtype *J, const realtype *dxdotdp);
 extern void xBdot_model_dirac(realtype *xBdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *xB, const realtype *w, const realtype *dwdx);
 extern void xdot_model_dirac(realtype *xdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w);
-extern void y_model_dirac(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
+extern void y_model_dirac(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w);
 
 class Model_model_dirac : public amici::Model_ODE {
 public:
@@ -62,7 +62,7 @@ class Model_model_dirac : public amici::Model_ODE {
 
     virtual amici::Model* clone() const override { return new Model_model_dirac(*this); };
 
-    const  char* getAmiciVersion() const { return "acb0375e4e0c93fb725d81b42e7e818c9ba7b65d"; };
+    const  char* getAmiciVersion() const { return ""; };
 
     virtual void fJ(realtype *J, const realtype t, const realtype *x, const double *p, const double *k, const realtype *h, const realtype *w, const realtype *dwdx) override {
         J_model_dirac(J, t, x, p, k, h, w, dwdx);
@@ -158,11 +158,11 @@ class Model_model_dirac : public amici::Model_ODE {
         dxdotdp_model_dirac(dxdotdp, t, x, p, k, h, ip, w, dwdp);
     }
 
-    virtual void fdydp(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) override {
+    virtual void fdydp(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp) override {
     }
 
-    virtual void fdydx(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {
-        dydx_model_dirac(dydx, t, x, p, k, h);
+    virtual void fdydx(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx) override {
+        dydx_model_dirac(dydx, t, x, p, k, h, w, dwdx);
     }
 
     virtual void fdzdp(double *dzdp, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) override {
@@ -220,8 +220,8 @@ class Model_model_dirac : public amici::Model_ODE {
         xdot_model_dirac(xdot, t, x, p, k, h, w);
     }
 
-    virtual void fy(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {
-        y_model_dirac(y, t, x, p, k, h);
+    virtual void fy(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w) override {
+        y_model_dirac(y, t, x, p, k, h, w);
     }
 
     virtual void fz(double *z, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {

diff --git a/models/model_dirac/model_dirac_dydx.cpp b/models/model_dirac/model_dirac_dydx.cpp
@@ -6,7 +6,7 @@ typedef amici::realtype realtype;
 
 using namespace amici;
 
-void dydx_model_dirac(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+void dydx_model_dirac(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx) {
   dydx[0+1*1] = 1.0;
 }
 
diff --git a/models/model_dirac/model_dirac_y.cpp b/models/model_dirac/model_dirac_y.cpp
@@ -6,7 +6,7 @@ typedef amici::realtype realtype;
 
 using namespace amici;
 
-void y_model_dirac(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+void y_model_dirac(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w) {
   y[0] = x[1];
 }
 
diff --git a/models/model_events/model_events.h b/models/model_events/model_events.h
@@ -1,6 +1,6 @@
 #ifndef _amici_model_events_h
 #define _amici_model_events_h
-/* Generated by amiwrap (R2017b) acb0375e4e0c93fb725d81b42e7e818c9ba7b65d */
+/* Generated by amiwrap (R2017b)  */
 #include <cmath>
 #include <memory>
 #include "amici/defines.h"
@@ -32,8 +32,8 @@ extern void dJzdz_model_events(double *dJzdz, const int iz, const realtype *p, c
 extern void deltasx_model_events(double *deltasx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const int ip, const int ie, const realtype *xdot, const realtype *xdot_old, const realtype *sx, const realtype *stau);
 extern void drzdx_model_events(double *drzdx, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
 extern void dxdotdp_model_events(realtype *dxdotdp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp);
-extern void dydp_model_events(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip);
-extern void dydx_model_events(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
+extern void dydp_model_events(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp);
+extern void dydx_model_events(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx);
 extern void dzdx_model_events(double *dzdx, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
 extern void qBdot_model_events(realtype *qBdot, const int ip, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *xB, const realtype *w, const realtype *dwdp);
 extern void root_model_events(realtype *root, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
@@ -47,7 +47,7 @@ extern void sz_model_events(double *sz, const int ie, const realtype t, const re
 extern void x0_model_events(realtype *x0, const realtype t, const realtype *p, const realtype *k);
 extern void xBdot_model_events(realtype *xBdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *xB, const realtype *w, const realtype *dwdx);
 extern void xdot_model_events(realtype *xdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w);
-extern void y_model_events(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
+extern void y_model_events(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w);
 extern void z_model_events(double *z, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h);
 
 class Model_model_events : public amici::Model_ODE {
@@ -76,7 +76,7 @@ class Model_model_events : public amici::Model_ODE {
 
     virtual amici::Model* clone() const override { return new Model_model_events(*this); };
 
-    const  char* getAmiciVersion() const { return "acb0375e4e0c93fb725d81b42e7e818c9ba7b65d"; };
+    const  char* getAmiciVersion() const { return ""; };
 
     virtual void fJ(realtype *J, const realtype t, const realtype *x, const double *p, const double *k, const realtype *h, const realtype *w, const realtype *dwdx) override {
         J_model_events(J, t, x, p, k, h, w, dwdx);
@@ -178,12 +178,12 @@ class Model_model_events : public amici::Model_ODE {
         dxdotdp_model_events(dxdotdp, t, x, p, k, h, ip, w, dwdp);
     }
 
-    virtual void fdydp(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) override {
-        dydp_model_events(dydp, t, x, p, k, h, ip);
+    virtual void fdydp(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp) override {
+        dydp_model_events(dydp, t, x, p, k, h, ip, w, dwdp);
     }
 
-    virtual void fdydx(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {
-        dydx_model_events(dydx, t, x, p, k, h);
+    virtual void fdydx(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx) override {
+        dydx_model_events(dydx, t, x, p, k, h, w, dwdx);
     }
 
     virtual void fdzdp(double *dzdp, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) override {
@@ -247,8 +247,8 @@ class Model_model_events : public amici::Model_ODE {
         xdot_model_events(xdot, t, x, p, k, h, w);
     }
 
-    virtual void fy(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {
-        y_model_events(y, t, x, p, k, h);
+    virtual void fy(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w) override {
+        y_model_events(y, t, x, p, k, h, w);
     }
 
     virtual void fz(double *z, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) override {

diff --git a/models/model_events/model_events_dydp.cpp b/models/model_events/model_events_dydp.cpp
@@ -6,7 +6,7 @@ typedef amici::realtype realtype;
 
 using namespace amici;
 
-void dydp_model_events(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip) {
+void dydp_model_events(double *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const int ip, const realtype *w, const realtype *dwdp) {
 switch (ip) {
   case 3: {
   dydp[0] = x[0]+x[1]+x[2];

diff --git a/models/model_events/model_events_dydx.cpp b/models/model_events/model_events_dydx.cpp
@@ -6,7 +6,7 @@ typedef amici::realtype realtype;
 
 using namespace amici;
 
-void dydx_model_events(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+void dydx_model_events(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx) {
   dydx[0+0*1] = p[3];
   dydx[0+1*1] = p[3];
   dydx[0+2*1] = p[3];

diff --git a/models/model_events/model_events_y.cpp b/models/model_events/model_events_y.cpp
@@ -6,7 +6,7 @@ typedef amici::realtype realtype;
 
 using namespace amici;
 
-void y_model_events(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h) {
+void y_model_events(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w) {
   y[0] = p[3]*(x[0]+x[1]+x[2]);
 }