kinetic2.c

GESPI 2.0动态系统模拟工具

#include "copyleft.h"

/*
    GEPASI - a simulator of metabolic pathways and other dynamical systems
    Copyright (C) 1989, 1992  Pedro Mendes
*/

/*************************************/
/*                                   */
/*      reaction rate equations      */
/*                                   */
/*          MICROSOFT C 6.00         */
/*           QuickC/WIN 1.0          */
/*             ULTRIX cc             */
/*              GNU gcc              */
/*                                   */
/*   (include here compilers that    */
/*   compiled GEPASI successfully)   */
/*                                   */
/*************************************/

#include <math.h>
#include "globals.h"
#include "globvar.h"
#include "datab.h"

/*****************************************/
/*                                       */
/*              CONVENTIONS              */
/*                                       */
/* all functions take as arguments:      */
/* double s[] - vector of concentrations */
/* int r - the index of the reaction     */
/*                                       */
/*     all fuctions return a double.     */
/*                                       */
/*                                       */
/*****************************************/


/*                               */
/* Alosteric Inhibition          */
/* As in Hofmeyr & van der Merwe */
/* CABIOS vol 2 243-249          */
/*                               */

double alinh( double s[], int r )
          /* eff[r][0] - substrate           */
          /* eff[r][1] - product             */
          /* eff[r][2] - inhibitor           */
          /* params[r][0] - Km for substrate    */
          /* params[r][1] - Km for product      */
          /* params[r][2] - max forward vel.    */
          /* params[r][3] - max reverse vel.    */
          /* params[r][4] - inhibition constant */
          /* params[r][5] - hill coeff (inhib)  */
{
 double nsb, npr;

 nsb = s[eff[r][0]]/params[r][0];
 npr = s[eff[r][1]]/params[r][1];
 return (nsb*params[r][2] - npr*params[r][3]) / (1 + nsb + npr + pow(s[eff[r][2]], params[r][5])/params[r][4]);
}



/*                               */
/* UNI-BI SEQUENTIAL             */
/* ( A = P + Q )                 */
/* functions deduced with reduce */
/* file A-PQ.REQ                 */
/*                               */
double unibisa( double s[], int r )
          /* eff[r][0] - A                   */
          /* eff[r][1] - P                   */
          /* eff[r][2] - Q                   */
          /* params[r][0] - KA                  */
          /* params[r][1] - KP                  */
          /* params[r][2] - KQ                  */
          /* params[r][3] - KPQ                 */
          /* params[r][4] - max forward vel.    */
          /* params[r][5] - max reverse vel.    */
{
 double na, np, npq;

 na  = s[eff[r][0]]/params[r][0];
 np  = s[eff[r][1]]/params[r][1];
 npq = s[eff[r][1]]*s[eff[r][2]]/params[r][3];

 return (na*params[r][4] - npq*params[r][5]) /
        ( 1 +
          na +
          na*np*(1+params[r][4]*(1-params[r][1]*params[r][2]/params[r][3])/params[r][5]) +
          npq +
          np +
          s[eff[r][2]]/params[r][2]
        );
}

double dunibisaa( double s[], int r, int e )
{
 double kpq2, kpkq, kpkpq, vrkpq, c;

 kpq2  = params[r][3]*params[r][3];
 kpkq  = params[r][1]*params[r][2];
 kpkpq = params[r][1]*params[r][3];
 vrkpq = params[r][3]*params[r][5];
 c     = (((kpkq-params[r][3])*params[r][4]-vrkpq)*params[r][2]*s[eff[r][1]]-params[r][5]*kpq2)*s[eff[r][0]]
         -(s[eff[r][2]]+params[r][2])*s[eff[r][1]]*params[r][0]*vrkpq
         -s[eff[r][2]]*params[r][0]*params[r][5]*kpkpq
         -params[r][0]*params[r][5]*kpq2;

 return (s[eff[r][1]]*params[r][2]+params[r][3])*params[r][0]*params[r][5]*params[r][5]*params[r][3]*
        (s[eff[r][1]]*s[eff[r][2]]*(params[r][3]*(params[r][4]+params[r][5])+params[r][4]*kpkq)+
         s[eff[r][2]]*params[r][4]*kpkpq + kpq2*params[r][5] )
        /
        c*c ;
}

double dunibisap( double s[], int r, int e )
{
 double kpq2, kpkq, kpkpq, vrkpq, c;

 kpq2  = params[r][3]*params[r][3];
 kpkq  = params[r][1]*params[r][2];
 kpkpq = params[r][1]*params[r][3];
 vrkpq = params[r][3]*params[r][5];
 c     = (((kpkq-params[r][3])*params[r][4]-vrkpq)*params[r][2]*s[eff[r][1]]-params[r][5]*kpq2)*s[eff[r][0]]
         -(s[eff[r][2]]+params[r][2])*s[eff[r][1]]*params[r][0]*vrkpq
         -s[eff[r][2]]*params[r][0]*params[r][5]*kpkpq
         -params[r][0]*params[r][5]*kpq2;

 return - (s[eff[r][0]]*params[r][2]*params[r][4]+s[eff[r][2]]*params[r][0]*params[r][5])*params[r][5]*kpq2*
          (s[eff[r][0]]*(params[r][3]*(params[r][4]+params[r][5])+params[r][4]*kpkq)+
           s[eff[r][2]]*params[r][0]*params[r][5]*params[r][1] + params[r][0]*vrkpq )
          /
          c*c ;
}


double dunibisaq( double s[], int r, int e )
{
 double kpq2, kpkq, kpkpq, vrkpq, c;

 kpq2  = params[r][3]*params[r][3];
 kpkq  = params[r][1]*params[r][2];
 kpkpq = params[r][1]*params[r][3];
 vrkpq = params[r][3]*params[r][5];
 c     = (((kpkq-params[r][3])*params[r][4]-vrkpq)*params[r][2]*s[eff[r][1]]-params[r][5]*kpq2)*s[eff[r][0]]
         -(s[eff[r][2]]+params[r][2])*s[eff[r][1]]*params[r][0]*vrkpq
         -s[eff[r][2]]*params[r][0]*params[r][5]*kpkpq
         -params[r][0]*params[r][5]*kpq2;

 return - (s[eff[r][1]]*params[r][2]+params[r][3])*params[r][0]*params[r][5]*params[r][5]*params[r][3]*
          (s[eff[r][0]]*s[eff[r][1]]*(params[r][3]*(params[r][4]+params[r][5])+params[r][4]*kpkq)+
           s[eff[r][0]]*params[r][4]*kpkpq + s[eff[r][1]]*params[r][0]*vrkpq )
          /
          c*c ;
}


/*                               */
/* UNI-BI                        */
/* ( A = 2P )                    */
/* functions deduced with reduce */
/* file A-PP.REQ                 */
/*                               */
double unibisb( double s[], int r )
          /* eff[r][0] - A                   */
          /* eff[r][1] - P                   */
          /* params[r][0] - KA                  */
          /* params[r][1] - KP                  */
          /* params[r][2] - max forward vel.    */
          /* params[r][3] - max reverse vel.    */
{
 double p2, vf2, vr_pl_2vf;

 p2        = s[eff[r][1]]*s[eff[r][1]];
 vf2       = params[r][2]*params[r][2];
 vr_pl_2vf = params[r][3] + 2*params[r][2];

 return 4*vf2*(vr_pl_2vf*s[eff[r][0]]*params[r][1]*params[r][1] - 2*p2*params[r][3]*params[r][0])
        /
        ( params[r][1]*( (s[eff[r][1]]*params[r][3]+4*params[r][2])*vr_pl_2vf*s[eff[r][0]]*params[r][1]
                  + 4*vr_pl_2vf*params[r][2]*params[r][0]*params[r][1]
                  + 8*(params[r][3]+params[r][2])*s[eff[r][1]]*params[r][2]*params[r][0]
                 ) + 8*p2*vf2*params[r][0]
        );
}

double dunibisba( double s[], int r, int e )
{
 double p2, vf2, vr_pl_2vf, vf_pl_vr, kp2, p_pl_kp, c;

 p2        = s[eff[r][1]]*s[eff[r][1]];
 vf2       = params[r][2]*params[r][2];
 vr_pl_2vf = params[r][3] + 2*params[r][2];
 vf_pl_vr  = params[r][2] + params[r][3];
 kp2       = params[r][1]*params[r][1];
 p_pl_kp   = s[eff[r][1]] + params[r][1];
 c         = params[r][1]*( (s[eff[r][1]]*params[r][3]+4*params[r][2])*vr_pl_2vf*s[eff[r][0]]*params[r][1]
                      + 4*vr_pl_2vf*params[r][2]*params[r][0]*params[r][1]
                      + 8*vf_pl_vr*s[eff[r][1]]*params[r][2]*params[r][0]
                   ) + 8*p2*vf2*params[r][0];

 return 8*(p2*s[eff[r][1]]*params[r][3]*params[r][3] + 2*params[r][2]*vf_pl_vr*(p2 + p_pl_kp*p_pl_kp))
        *vr_pl_2vf*vf2*params[r][0]*kp2
        / c*c;
}

double dunibisbp( double s[], int r, int e )
{
 double p2, vf2, vr2, vr_pl_2vf, vf_pl_vr, kakp, vrvf, vfvr_big, c;

 p2        = s[eff[r][1]]*s[eff[r][1]];
 vf2       = params[r][2]*params[r][2];
 vr2       = params[r][3]*params[r][3];
 vr_pl_2vf = params[r][3] + 2*params[r][2];
 vf_pl_vr  = params[r][2] + params[r][3];
 kakp      = params[r][0]*params[r][1];
 vrvf      = params[r][2]*params[r][3];
 vfvr_big  = vr2 + 3*vrvf + 2*vf2;
 c         = params[r][1]*( (s[eff[r][1]]*params[r][3]+4*params[r][2])*vr_pl_2vf*s[eff[r][0]]*params[r][1]
                      + 4*vr_pl_2vf*params[r][2]*kakp
                      + 8*vf_pl_vr*s[eff[r][1]]*params[r][2]*params[r][0]
                   ) + 8*p2*vf2*params[r][0];

 return -4*vf2*params[r][1]*( s[eff[r][0]]*s[eff[r][0]]*params[r][1]*params[r][1]*params[r][1]*params[r][3]*(vr2 + 4*params[r][2]*vf_pl_vr)
                      + 2*s[eff[r][0]]*p2*kakp*vr2*vr_pl_2vf
                      + 16*s[eff[r][0]]*s[eff[r][1]]*kakp*params[r][2]*vfvr_big
                      + 8*s[eff[r][0]]*kakp*params[r][1]*params[r][2]*vfvr_big
                      + 16*p2*params[r][0]*params[r][0]*vrvf*vf_pl_vr
                      + 16*s[eff[r][1]]*kakp*params[r][0]*vrvf*vf_pl_vr
                    )
        / c*c ;
}