/* 
 * This program is free software; you can redistribute it and/or modify 
 * it under the terms of the GNU General Public License as published by 
 * the Free Software Foundation; either version 2 of the License, or 
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful, but 
 * WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 
 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 
 * for more details.
 * 
 * You should have received a copy of the GNU General Public License along 
 * with this program; if not, write to the Free Software Foundation, Inc., 
 * 59
 */

#include <iostream>
#include <math.h>
# include "BiStableSwitch31st.h"

using namespace std;

            double kAssHSLLuxR			= 0.0001;
            double kDissHSLLuxR	                = 0.003;
            double kMaxProductionRNA            = 0.44;
            double kHalfMaxProductionRNA        = 0.013;
            double nplasmids                    = 10;

            double kDegRNA                      = 0.00288;
            double kDegLuxI			= 0.00288;
            double kDegLuxR	                = 0.001;
            double kDegHSL                      = 0.000187;
            double kDegGFP                      = 0.00385;

            double kConvHSL	                = 0.03 ;
            double n                            = 2;
            double kdiff                        = 0.4 ;
            double kTranslation                 = 0.1 ;
            double kDimer                       = 3 * pow(10,-5) ;
            double kDissHSLLuxR2                = 0.1 ;


double *ode(double y[], double dy[]) {



/// Section 0: HSLOut Concentration,

dy[0] =  - kDegHSL * y[0] - ( y[15] * ( y[11]/y[12] ) * ( y[0] - y[1] ) )  - ( y[0] * dy[10] / y[11] ) ;

/// Section 1: HSL Concentration.

dy[1] = ( y[14] * y[7] ) * ( y[11]/y[12] ) - kAssHSLLuxR * y[1] * y[7] + 2 * kDissHSLLuxR * y[2] - kDegHSL * y[1] + kDegLuxR  * y[2] + ( y[15] * ( y[0] - y[1] ) ) - ( y[1] * dy[10] / y[11] ) ;

/// Section 2: HSL-LuxR Concentration.

dy[2] = ( kAssHSLLuxR * y[1] * y[8] ) - ( kDissHSLLuxR * y[2] ) - ( kDegLuxR * y[2] ) - ( kDegHSL * y[2] ) - ( kDimer *y[2] ) + ( 2 * kDissHSLLuxR2 * y[3] ) - ( y[2] * dy[10] / y[11] );

/// Section 3: HSL-LuxR^2 Concentration.

dy[3] =  ( kDimer * y[2] ) - ( kDissHSLLuxR2 * y[3] ) - ( kDegHSL * y[3] ) - ( kDegLuxR * y[3] ) - ( y[3] * dy[10] / y[11] ) ;

/// Section 4: LuxI mRNA Concentration.

dy[4] = ( y[13] * ( kMaxProductionRNA * ( pow(y[3],n ))/( pow(kHalfMaxProductionRNA, n) + pow(y[3], n) ) ) ) - ( kDegRNA * y[4] )  - ( y[4] * dy[10] / y[11] ) ;

/// Section 5: LuxR mRNA Concentration.

 dy[5] = ( y[13] * ( kMaxProductionRNA * ( pow(y[3],n ))/( pow(kHalfMaxProductionRNA, n) + pow(y[3], n) ) ) )  - ( kDegRNA * y[5] )  - ( y[5] * dy[10] / y[11] ) ;

 /// Section 6: GFP mRNA Concentration.

 dy[6] = ( y[13] * ( kMaxProductionRNA * ( pow(y[3],n ))/( pow(kHalfMaxProductionRNA, n) + pow(y[3], n) ) ) ) - ( kDegRNA * y[6] )  - ( y[6] * dy[10] / y[11] ) ;

/// Section 7: LuxI Concentration.

dy[7] = ( kTranslation * y[4] ) - (kDegLuxI * y[7] )  - ( y[7] * dy[10] / y[11] ) ;

/// Section 8: LuxR Concentration.

dy[8] = ( kTranslation * y[5] ) - ( kAssHSLLuxR * y[1] * y[8] ) + ( kDissHSLLuxR * y[2] ) + ( kDegHSL * y[2] ) - ( kDegLuxR * y[8] ) - ( y[8] * dy[10] / y[11] ) ;

/// Section 9: GFP Concentration.

dy[9] = ( kTranslation * y[6] ) - ( kDegGFP * y[9] ) - ( y[9] * dy[10] / y[11] ) ;

return dy;

}