# with the following commands at the prompt in the OpenCOR Python console: # # In [1]: cd path/to/folder_this_file_is_in # In [2]: run Fig1_sim.py import opencor as oc import matplotlib.pyplot as plt import numpy as np # The prefix of the saved output file name V_var = 0; V_hold = 1; if V_var: prefilename = 'Fig5' # Load the simulation file simfile = 'Tong_2011.sedml' simulation = oc.open_simulation(simfile) data = simulation.data() # Reset states variables and parameters # simulation.reset(True) # Set constant values start = 0 end = 1000 pointInterval = 0.5 data.set_starting_point(start) data.set_ending_point(end) data.set_point_interval(pointInterval) # simulation.reset(True) data.constants()['interface/membrane_potential/v_init'] = -80 data.constants() ['interface/parameters/kmca'] = 0.001 data.constants() [ 'interface/I_h/a_correction' ] = 1 simulation.run() # Access simulation results results = simulation.results() # Data to save varName = np.array(["Time", "qss", "rss", "qtc", "r1tc","r2tc", "ik1", "v"]) vars = np.reshape(varName, (1,8)) rows = end * 2 + 1 print(rows) # Access simulation results results = simulation.results() # Grab a some algebraic variable results r = np.zeros((rows, len(varName))) r [ :, 0 ] = results.voi().values() r [ :, 1 ] = results.algebraic()['interface/I_K1/qss'].values()#q2 r [ :, 2 ] = results.algebraic()['interface/I_K1/rss'].values() r [ :, 3 ] = results.algebraic() [ 'interface/I_K1/qtc'].values() r [ :, 4 ] = results.algebraic() [ 'interface/I_K1/r1tc'].values() r [ :, 5 ] = results.algebraic() [ 'interface/I_K1/r2tc'].values() r [ :, 6 ] = results.algebraic() [ 'interface/I_K1/ik1'].values() r [ :, 7 ] = results.states() [ 'interface/membrane_potential/v' ].values() filename = '%s.csv' % (prefilename) np.savetxt(filename, vars, fmt='%s', delimiter=",") with open(filename, "ab") as f: np.savetxt(f, r, delimiter=",") f.close elif V_hold: prefilename = 'Fig5_4' # Load the simulation file simfile = 'VoltageVAR_Tong_2011.sedml' simulation = oc.open_simulation(simfile) data = simulation.data() simulation.reset(True) # Reset states variables and parameters #simulation.reset(True) # Set constant values start = 0 end = 10000 pointInterval = 1 data.set_starting_point(start) data.set_ending_point(end) data.set_point_interval(pointInterval) time_clamp = 40 v_clamp = [ -30, -20, -10, 0, 10, 20, 30, 40] v_hold = [-380, -40] simulation.reset(True) for j in range(len(v_hold)): data.constants()['interface/membrane_potential/v_hold'] = v_hold[j] for i in range(len(v_clamp)): data.constants() [ 'interface/membrane_potential/v_clamp' ] = v_clamp [ i ] print(i,j) data.constants() [ 'interface/membrane_potential/time_clamp' ] = time_clamp data.constants() [ 'interface/parameters/ki' ] = 110 data.constants() [ 'interface/parameters/ko' ] = 4 #data.constants() [ 'interface/membrane_potential/v_clamp' ] = v[ i ] simulation.run() # Access simulation results results = simulation.results() # Data to save # Data to save varName = np.array([ "Time", "ik1", "v","max" ]) vars = np.reshape(varName, (1, 4)) rows = end * 1 + 1 # Grab a some algebraic variable results r = np.zeros((rows, len(varName))) r [ :, 0 ] = results.voi().values() r [ :, 1 ] = results.algebraic() [ 'interface/I_K1/ik1' ].values() r [ :, 2 ] = results.algebraic() [ 'interface/membrane_potential/v' ].values() r[:,3] = max (r [ :, 1 ]) # Save the simulation result of the last run filename = '%s_%s_%s.csv' % (prefilename, j,i) np.savetxt(filename, vars, fmt='%s', delimiter=",") with open(filename, "ab") as f: np.savetxt(f, r, delimiter=",") f.close # simulation.reset(True) # # # simulation.run() # # Access simulation results # results = simulation.results() # # Data to save # varName = np.array([ "Time", "ik1", "v" ]) # vars = np.reshape(varName, (1, 3)) # rows = end * 2 + 1 # # # Grab a some algebraic variable results # r = np.zeros((rows, len(varName))) # r [ :, 0 ] = results.voi().values() # r [ :, 1 ] = results.algebraic() [ 'interface/I_K1/ik1'].values() # r [ :, 2 ] = results.algebraic() [ 'interface/membrane_potential/v' ].values() # # # Save the simulation result of the last run # filename = '%s_%s.csv' % (prefilename, i) # np.savetxt(filename, vars, fmt='%s', delimiter=",") # with open(filename, "ab") as f: # np.savetxt(f, r, delimiter=",") # f.close else: prefilename = 'Fig5_5' # Load the simulation file simfile = 'VoltageVAR_Tong_2011.sedml' simulation = oc.open_simulation(simfile) data = simulation.data() # Reset states variables and parameters simulation.reset(True) # Set constant values start = 0 end = 10000 pointInterval = 0.5 data.set_starting_point(start) data.set_ending_point(end) data.set_point_interval(pointInterval) time_clamp = 40 v_clamp = [-40,-30,-20,-10,0, 10 ] for i in range(len(v_clamp)): data.constants() [ 'interface/membrane_potential/v_hold' ] = -80 simulation.reset(True) data.constants() [ 'interface/parameters/kmca' ] = 0.001 data.constants() [ 'interface/membrane_potential/v_clamp' ] = v_clamp [ i ] data.constants() [ 'interface/membrane_potential/time_clamp' ] = time_clamp simulation.run() # Access simulation results results = simulation.results() # Data to save # Data to save varName = np.array([ "Time", "ik1", "v" ]) vars = np.reshape(varName, (1, 3)) rows = end * 2 + 1 # Grab a some algebraic variable results r = np.zeros((rows, len(varName))) r [ :, 0 ] = results.voi().values() r [ :, 1 ] = results.algebraic() ['interface/I_K1/ik1'].values() r [ :, 2 ] = results.algebraic() [ 'interface/membrane_potential/v'].values() # Save the simulation result of the last run filename = '%s_%s.csv' % (prefilename, i) np.savetxt(filename, vars, fmt='%s', delimiter=",") with open(filename, "ab") as f: np.savetxt(f, r, delimiter=",") f.close