# To reproduce Figure 1 in the associated Physiome paper,
# execute this script from the command line:
#
#   cd [PathToThisFile]
#   [PathToOpenCOR]/pythonshell Figure7.py

import matplotlib.pyplot as plt
import opencor as opencor
import numpy as np

simulation = opencor.open_simulation("model.sedml")
data = simulation.data()
data.set_ending_point(3000)
data.set_point_interval(10)




def run_sim1(gl_l):
    simulation.reset(True)
    simulation.clear_results()
    data.constants()["Cell_concentration/L_A"] = 6e-5
    data.constants()["Cell_concentration/L_B"] = 6e-5
    data.constants()["Blood_concentrations/v_B"] = 1e-17
    data.constants()["Blood_concentrations/glucose_in"] = 0.004
    data.constants()["A_GLUT2/n_GLUT"] = 1e8
    data.constants()["Blood_concentrations/Q_in"] = 1e-17
    # data.constants()["Blood_concentrations/v_w1"] = 1.8e-4
    data.constants()["phenomonological_constants/n_SGLT"] = 3e7
    data.constants()["parameters/k0_12"] = 12000
    data.constants()["parameters/k0_61"] = 15
    data.constants()["Apical_concentrations/glucose_m"] = gl_l
    simulation.run()
    ds = simulation.results().data_store()
    v_cell = ds.voi_and_variables()["Cell_concentration/v_cell"].values()[-1]

    return (v_cell)


def run_sim2(gl_l):
    simulation.reset(True)
    simulation.clear_results()
    data.constants()["Cell_concentration/L_A"] = 6e-5
    data.constants()["Cell_concentration/L_B"] = 6e-5
    data.constants()["Blood_concentrations/v_B"] = 1e-17
    data.constants()["Blood_concentrations/glucose_in"] = 0.004
    data.constants()["A_GLUT2/n_GLUT"] = 0
    data.constants()["Blood_concentrations/Q_in"] = 1e-17
    # data.constants()["Blood_concentrations/v_w1"] = 1.8e-4
    data.constants()["phenomonological_constants/n_SGLT"] = 3e7
    data.constants()["parameters/k0_12"] = 12000
    data.constants()["parameters/k0_61"] = 15
    data.constants()["Apical_concentrations/glucose_m"] = gl_l
    simulation.run()
    ds = simulation.results().data_store()
    v_cell = ds.voi_and_variables()["Cell_concentration/v_cell"].values()[-1]

    return (v_cell)





#
#
if __name__ == '__main__':

    #
    X = [10,15,20,25,30,35,40,45,50]
    gl_l = [0.01, 0.015, 0.02, 0.025, 0.03, 0.035, 0.04, 0.045, 0.05]

    y_label = ["Blood Glucose", "Cell Glucose", "GLUT2 Flux"]

    plt.figure(figsize=(10, 6))
    plt.subplot(1, 2, 1)
    Gl_l_results = np.empty(len(gl_l))
    for i in range(len(gl_l)):
        Gl_l = run_sim1(gl_l[i])
        Gl_l_results[i] = Gl_l
    print(Gl_l_results)
    results = np.empty(len(Gl_l_results))
    for i in range(len(Gl_l_results)-1):
        result = 1-((Gl_l_results[0]-Gl_l_results[i+1])/Gl_l_results[0])
        results[i] = result
    print(results)

    final_results = np.array(1)
    for i in range(len(results) - 1):
        final_results = np.append(final_results, results[i])


    plt.plot(X, final_results, 'blue', linewidth=2.5, marker='d', markevery=slice(0, 5000, 100),
             label='With Apical GLUT2')

    Gl_l_results = np.empty(len(gl_l))
    for i in range(len(gl_l)):
        Gl_l = run_sim2(gl_l[i])
        Gl_l_results[i] = Gl_l
    print(Gl_l_results)
    results = np.empty(len(Gl_l_results))
    for i in range(len(Gl_l_results) - 1):
        result = 1 - ((Gl_l_results[0] - Gl_l_results[i + 1]) / Gl_l_results[0])
        results[i] = result
    print(results)

    final_results = np.array(1)
    for i in range(len(results) - 1):
        final_results = np.append(final_results, results[i])

    plt.plot(X, final_results, 'red', linewidth=2.5, marker='d', markevery=slice(0, 5000, 100),
             label='Without Apical GLUT2')

    plt.grid()
    plt.xlim(10, 50)
    plt.tick_params(axis='x', which='major', labelsize=10)
    plt.tick_params(axis='y', which='major', labelsize=10)
    # ~ plt.ylim(0.0, 0.2)
    plt.xlabel('Luminal Glucose (mM)', fontsize='10')
    plt.ylabel('Cell Volume (%)', fontsize='10')
    plt.legend(loc='best', fontsize='10')
    plt.title('A')

    plt.subplot(1, 2, 2)

    v_cell_Naftalin_with = np.array((7.56E-17, 8E-17, 8.40E-17, 8.65E-17, 9.07E-17))
    v_cell_Naftalin_without = np.array((8.40E-17, 9.72E-17, 1.12E-16, 1.29E-16, 1.56E-16))

    results = np.empty(len(v_cell_Naftalin_with))
    for i in range(len(v_cell_Naftalin_with) - 1):
        result = 1 - ((v_cell_Naftalin_with[0] - v_cell_Naftalin_with[i+1]) / v_cell_Naftalin_with[0])
        results[i] = result
    print(results)
    final_results = np.array(1)
    for i in range(len(results) - 1):
        final_results = np.append(final_results, results[i])
    X1 = [10, 20, 30, 40, 50]
    plt.plot(X1, final_results, 'blue', linewidth=2.5, marker='d', markevery=slice(0, 5000, 100),
             label='With Apical GLUT2')

    results = np.empty(len(v_cell_Naftalin_without))
    for i in range(len(v_cell_Naftalin_without) - 1):
        result = 1 - ((v_cell_Naftalin_without[0] - v_cell_Naftalin_without[i + 1]) / v_cell_Naftalin_without[0])
        results[i] = result
    print(results)
    final_results = np.array(1)
    for i in range(len(results) - 1):
        final_results = np.append(final_results, results[i])

    plt.plot(X1, final_results, 'red', linewidth=2.5, marker='d', markevery=slice(0, 5000, 100),
             label='Without Apical GLUT2')

    plt.grid()
    plt.xlim(10, 50)
    plt.tick_params(axis='x', which='major', labelsize=10)
    plt.tick_params(axis='y', which='major', labelsize=10)
    # ~ plt.ylim(0.0, 0.2)
    plt.xlabel('Luminal Glucose (mM)', fontsize='10')
    plt.ylabel('Cell Volume (%)', fontsize='10')
    plt.legend(loc='best', fontsize='10')
    plt.title('B')


    plt.savefig('Figure07.png')
    plt.show()