def test_truth_table(testmodel, input_ids, output_ids, truth_table,
         ht=0.8, lt=0.2, delay=99, plot=False):
    import tellurium as te
    r = te.loada(testmodel)

    sims = []

    for row in truth_table:
        message = ['When']
        for i, input_val in enumerate(row[0]):
            message.append(input_ids[i] + ':' + str(input_val))
            r[input_ids[i]] = input_val

        sim = r.simulate(0, delay + 1, delay + 1, ['time'] + input_ids + output_ids)
        sims.append(sim)

        for i, output_val in enumerate(row[1]):
            offset = len(input_ids) + 1 # Time + length of inputs
            ind = i + offset
            full_message = ' '.join(message + [
                'then',
                output_ids[i] + ':' + str(output_val),
                '; Found %s = %f' % (output_ids[i], sim[delay][ind])
            ])
            print full_message
            if output_val == 0:
                assert sim[delay][ind] < lt, full_message
            else:
                assert sim[delay][ind] > ht, full_message
    return r, sims

    def run(self,func=None):
        """Allows the user to set the data from a File
        This data is to be compared with the simulated data in the process of parameter estimation
        
        Args:
            func: An Optional Variable with default value (None) which by default run differential evolution
                which is from scipy function. Users can provide reference to their defined function as argument.
            

        Returns:
            The Value of the parameter(s) which are estimated by the function provided.
        
        .. sectionauthor:: Shaik Asifullah <[email protected]>
        
        
        """
        
        self._parameter_names = self.bounds.keys()
        self._parameter_bounds = self.bounds.values()
        self._model_roadrunner = te.loada(self.model.model)
        x_data = self.data[:,0]
        y_data = self.data[:,1:]
        arguments = (x_data,y_data)

        if(func is not None):
            result = differential_evolution(self._SSE, self._parameter_bounds, args=arguments)
            return(result.x)
        else:
            result = func(self._SSE,self._parameter_bounds,args=arguments)
            return(result.x)

    def test_plot(self):
        """ Regression tests for plotting.
        The following calls should work. """
        r = te.loada("""
            S1 -> S2; k1*S1;
            k1 = 0.1; S1 = 40; S2 = 0.0;
        """)
        print(type(r))

        s = r.simulate(0, 100, 21)
        # no argument version
        r.plot()
        # plot with data
        r.plot(s)
        # plot with named data
        r.plot(result=s)
        # plot without legend
        r.plot(s)
        # plot without showing
        r.plot(s, show=False)
        r.plot(s, show=True)  # no show
        # plot with label, title, axis and legend
        r.plot(s, xlabel="x", ylabel="y", xlim=[0, 10], ylim=[0, 10], grid=True)
        # plot with additional plot settings from matplotlib
        r.plot(s, alpha=0.1, color="blue", linestyle="-", marker="o")

 def test_setSeed(self):
     r = te.loada("""
         S1 -> S2; k1*S1;
         k1 = 0.1; S1 = 40; S2 = 0.0;
     """)
     r.setSeed(123)
     self.assertEqual(123, r.getSeed())

 def test_getSeed(self):
     r = te.loada("""
         S1 -> S2; k1*S1;
         k1 = 0.1; S1 = 40; S2 = 0.0;
     """)
     seed = r.getSeed()
     self.assertIsNotNone(seed)

 def spark_work(model_with_parameters):
     import tellurium as te
     if(antimony == "antimony"):
         model_roadrunner = te.loada(model_with_parameters[0])
     else:
         model_roadrunner = te.loadSBMLModel(model_with_parameters[0])
     parameter_scan_initilisation = te.ParameterScan(model_roadrunner,**model_with_parameters[1])
     simulator = getattr(parameter_scan_initilisation, function_name)
     return(simulator())

 def test_loada(self):
     rr = te.loada('''
         model example0
           S1 -> S2; k1*S1
           S1 = 10
           S2 = 0
           k1 = 0.1
         end
     ''')
     self.assertIsNotNone(rr.getModel())

 def test_draw(self):
     r = te.loada("""
         S1 -> S2; k1*S1;
         k1 = 0.1; S1 = 40; S2 = 0.0;
     """)
     try:
         import pygraphviz
         r.draw()
     except ImportError:
         pass

 def fromAntimony(cls, antimonyStr, location, master=None):
     """ Create SBMLAsset from antimonyStr
     :param antimonyStr:
     :type antimonyStr:
     :param location:
     :type location:
     :return:
     :rtype:
     """
     r = te.loada(antimonyStr)
     raw = r.getSBML()
     return cls.fromRaw(raw=raw, location=location, filetype='sbml', master=master)

 def test_README_example(self):
     """ Tests the source example in the main README.md. """
     import tellurium as te
     rr = te.loada('''
         model example0
           S1 -> S2; k1*S1
           S1 = 10
           S2 = 0
           k1 = 0.1
         end
     ''')
     result = rr.simulate(0, 40, 500)
     te.plotArray(result)

 def fromAntimony(cls, antimonyStr, location, master=None):
     """ Create SBMLAsset from antimonyStr
     :param antimonyStr:
     :type antimonyStr:
     :param location:
     :type location:
     :return:
     :rtype:
     """
     warnings.warn('Use inline_omex instead.', DeprecationWarning)
     r = te.loada(antimonyStr)
     raw = r.getSBML()
     return cls.fromRaw(raw=raw, location=location, filetype='sbml', master=master)

 def stochastic_work(model_object):
     import tellurium as te
     if model_type == "antimony":
         model_roadrunner = te.loada(model_object.model)
     else:
         model_roadrunner = te.loadSBMLModel(model_object.model)
     model_roadrunner.integrator = model_object.integrator
     # seed the randint method with the current time
     random.seed()
     # it is now safe to use random.randint
     model_roadrunner.setSeed(random.randint(1000, 99999))
     model_roadrunner.integrator.variable_step_size = model_object.variable_step_size
     model_roadrunner.reset()
     simulated_data = model_roadrunner.simulate(model_object.from_time, model_object.to_time, model_object.step_points)
     return([simulated_data.colnames,np.array(simulated_data)])

    def test_seed(self):
        r = te.loada('''
        S1 -> S2; k1*S1; k1 = 0.1; S1 = 40; S2 = 0;
        ''')

        # Simulate from time zero to 40 time units
        result = r.gillespie(0, 40, 11)

        # Simulate on a grid with 10 points from start 0 to end time 40
        result = r.gillespie(0, 40, 10)

        # Simulate from time zero to 40 time units using the given selection list
        # This means that the first column will be time and the second column species S1
        result = r.gillespie(0, 40, 11, ['time', 'S1'])

        # Simulate from time zero to 40 time units, on a grid with 20 points
        # using the give selection list
        result = r.gillespie(0, 40, 20, ['time', 'S1'])

    def test_plot2DParameterScan(self):
        """Test plot2DParameterScan."""
        import tellurium as te
        from tellurium.analysis.parameterscan import plot2DParameterScan
        r = te.loada("""
        model test
           J0: S1 -> S2; Vmax * (S1/(Km+S1))
            S1 = 10; S2 = 0;
            Vmax = 1; Km = 0.5;
        end
        """)
        s = r.simulate(0, 50, 101)
        # r.plot(s)

        import numpy as np
        plot2DParameterScan(r,
                            p1='Vmax', p1Range=np.linspace(1, 10, num=5),
                            p2='Vmax', p2Range=np.linspace(0.1, 1.0, num=5),
                            start=0, end=50, points=101)

    def test_complex_simulation(self):
        """ Test complex simulation. """
        model = '''
        model feedback()
           // Reactions:
           J0: $X0 -> S1; (VM1 * (X0 - S1/Keq1))/(1 + X0 + S1 + S4^h);
           J1: S1 -> S2; (10 * S1 - 2 * S2) / (1 + S1 + S2);
           J2: S2 -> S3; (10 * S2 - 2 * S3) / (1 + S2 + S3);
           J3: S3 -> S4; (10 * S3 - 2 * S4) / (1 + S3 + S4);
           J4: S4 -> $X1; (V4 * S4) / (KS4 + S4);

          // Species initializations:
          S1 = 0; S2 = 0; S3 = 0;
          S4 = 0; X0 = 10; X1 = 0;

          // Variable initialization:
          VM1 = 10; Keq1 = 10; h = 10; V4 = 2.5; KS4 = 0.5;
        end
        '''
        r = te.loada(model)
        result = r.simulate(0, 40, 101)
        r.plot(result)

def simulate(antModelStr, constants):
    rr = te.loada(antModelStr.format(**constants))
    rr.simulate(0, 800000, 1000)
    rr.plot()
    rr.getSteadyStateValues()
    return rr.P

r = te.loada("""
    model normalized_species()
    
    # conversion between active (SA) and inactive (SI)
    J1: SA -> SI; k1*SA - k2*SI;
    k1 = 0.1; k2 = 0.02;
    
    # species
    species SA, SI, ST;
    SA = 10.0; SI = 0.0;
    const ST := SA + SI;
    
    SA is "active state S";
    SI is "inactive state S";
    ST is "total state S";
    
    # normalized species calculated via assignment rules
    species SA_f, SI_f;
    SA_f := SA/ST;
    SI_f := SI/ST;
    
    SA_f is "normalized active state S";
    SI_f is "normalized inactive state S";
    
    # parameters for your function
    P = 0.1;
    tau = 10.0;
    nA = 1.0;
    nI = 2.0;
    kA = 0.1;
    kI = 0.2;
    # now just use the normalized species in some math
    F := ( (1-(SI_f^nI)/(kI^nI+SI_f^nI)*(kI^nI+1) ) * ( (SA_f^nA)/(kA^nA+SA_f^nA)*(kA^nA+1) ) -P)*tau;

    end
""")

#Output requested: Amount of people in each state for years 1-10.
#
#This implementation uses a flattened model separated to M=Male and F=Female

import tellurium as te

r = te.loada ('''
MJ00: MH -> MS; MH*0.2
MJ02: MH -> MD; MH*0.01
MJ10: MS -> MH; MS*0.1
MJ12: MS -> MD; MS*0.3

FJ00: FH -> FS; FH*0.1
FJ02: FH -> FD; FH*0.01
FJ10: FS -> FH; FS*0.1
FJ12: FS -> FD; FS*0.3

MH = 50
MS = 0
MD = 0

FH = 50
FS = 0
FD = 0

''')

# This will create the SBML XML file
te.saveToFile ('Example3.xml', r.getSBML())

r.setIntegrator('gillespie')
r.integrator.variable_step_size = True

 k10 = random.uniform(0.0001, 1)
 k11 = random.uniform(0.0001, 1)
 k12 = random.uniform(1.0 / 2500.0, 1.0 / 10000000.0)
 Km9 = random.uniform(0.0001, 0.01)
 Km10 = random.uniform(0.0001, 0.01)
 Km11 = random.uniform(0.0001, 0.01)
 inter = 0.0
 S1 = 10.0    
 
 constants_ref = {'k1': str(k1), 'k2': str(k2), 'k3': str(k3), 'k4': str(k4), 
     'k5': str(k5), 'k6': str(k6), 'k7': str(k7), 'k9': str(k9),
     'k10': str(k10), 'k11': str(k11), 'k12': str(k12), 'Km9' : str(Km9), 
     'Km10': str(Km10), 'Km11': str(Km11), 'inter': str(inter), 'k8': str(k8),
     'S1': str(S1)}
     
 rr = te.loada(model.format(**constants))
 rr.getSteadyStateValues()
 S4_ref = rr.S4
 print 'S4_ref =', S4_ref
     
 inter = random.uniform(0.0000000001, 1.0)
 constants = {'k1': str(k1), 'k2': str(k2), 'k3': str(k3), 'k4': str(k4), 
     'k5': str(k5), 'k6': str(k6), 'k7': str(k7), 'k9': str(k9),
     'k10': str(k10), 'k11': str(k11), 'k12': str(k12), 'Km9' : str(Km9), 
     'Km10': str(Km10), 'Km11': str(Km11), 'inter': str(inter), 'k8': str(k8),
     'S1': str(S1)}
 
 
 rr = te.loada(model.format(**constants))
 rr.getSteadyStateValues()
 S4 = rr.S4

# -*- coding: utf-8 -*-
"""
Simple UniUni reaction with first-order mass-action kinetics.
"""
from __future__ import print_function
import tellurium as te

model = '''
model pathway()
    S1 -> S2; k1*S1

    # Initialize values
    S1 = 10; S2 = 0
    k1 = 1
end
'''

# load model
r = te.loada(model)

# carry out a time course simulation
# arguments are: start time, end time, number of points
result = r.simulate(0, 10, 100)

# plot results
r.plotWithLegend(result)