# -*- coding: utf-8 -*- """ Created on Fri May 10 15:13:35 2019 @author: ugupta """ import os os.environ["PATH"] += os.pathsep + 'C:/Program Files (x86)/Graphviz2.38/bin/' import pydot import numpy as np import pandas as pd from graphviz import * from subprocess import check_call from Reactions.reactions_import import * from Legend.legend import * # Inputs from user include species_comp, reactions file # along with cutoff rate, elements desired, rank sep and node sep if needed ''' Demonstration script for reading in reaction rates and reactions to generate visualizations ''' def get_color(surf_cov, phase): # print(surf_cov, phase) if phase == 'Gas': color = 'none' else: if float(surf_cov) >= 0.85: color = 'coral' elif 0.7 <= float(surf_cov) < 0.85: color = 'orangered' elif 0.55 <= float(surf_cov) < 0.7: color = 'red' elif 0.4 <= float(surf_cov) < 0.55: color = 'darkorange' elif 0.2 <= float(surf_cov) < 0.4: color = 'orange' elif 0.1 <= float(surf_cov) < 0.2: color = 'gold' elif 0.01 <= float(surf_cov) < 0.1: color = 'lawngreen' elif 1E-5 <= float(surf_cov) < 0.01: color = 'green' elif 1E-10 <= float(surf_cov) < 1E-5: color = 'dodgerblue' else: color = 'purple' return color def normalization_type(norm): if norm == 1: # print('Normalization requested is Max. Reaction Rate') return 'Normalization_MaxReactionRate.txt' elif norm == 2: # print('Normalization requested is Net Reaction Rate') return 'Normalization_NetReactionRate.txt' elif norm == 3: # print('Normalization requested is Local consumption') return 'Normalization_LocalConsumption.txt' else: print('No such normalization exists') def print_header(fname, Rank_Sep, Node_Sep): # print('Inside print header') # print(fname) f = open(fname, "w+") f.write('digraph G {\n') f.write('splines = true;\n') f.write('graph [bgcolor=lightgray, resolution=64, fontname=Arial, fontcolor=blue, fontsize=36];\n') f.write('node [fontsize=12];\n') f.write('edge [fontsize=30];\n') f.write('label = "Reaction Path Analysis";\n') f.write('labelloc = "t";\n') f.write('center=1;\n') f.write('size="10,10";\n') f.write('ranksep="' + str(Rank_Sep) + ' equally";\n') f.write('nodesep="' + str(Node_Sep) + ' equally";\n') f.write('rankdir=TB;\n') f.write('bgcolor=white;\n') f.close() def print_header_species(fname, Rank_Sep, Node_Sep): # print('Inside print header') # print(fname) f = open(fname, "w+") f.write('digraph G {\n') f.write('splines = true;\n') f.write('graph [bgcolor=lightgray, resolution=64, fontname=Arial, fontcolor=blue, fontsize=36];\n') f.write('node [fontsize=12];\n') f.write('edge [fontsize=30];\n') f.write('label = "Reaction Path Analysis";\n') f.write('labelloc = "t";\n') f.write('center=1;\n') f.write('size="10,10";\n') f.write('ranksep="' + str(Rank_Sep) + ' equally";\n') f.write('nodesep="' + str(Node_Sep) + ' equally";\n') f.write('rankdir=LR;\n') f.write('bgcolor=white;\n') f.close() def normalization(fname,normalization_rate): # print('IN here',openmkm_species_list) for i in range(len(openmkm_species_list)): if openmkm_species_consider[i] == True: count = 0 #Keeping a count of how many reactions with reaction rate more than cutoff are present for j in range(len(openmkm_reaction_strings)): #First checking if a reaction exists above the cutoff rate if float(stoich_matrix[i][j] != 0.0): if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): count = count + 1 #So if the species is being consumed in the reaction if float(rpa_local[i][j]) <= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: # print(openmkm_species_list[i],prod) #check if the product needs to be considered for visualization if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*20))) + int(1) arrowsize = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*120))) + int(1) edgelabel = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate))) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') else: linewidth = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*20))) + int(1) arrowsize = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*120))) + int(1) edgelabel = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate))) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') if float(rpa_local[i][j]) >= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: # print(openmkm_species_list[i], prod, 'new') if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*20))) + int(1) arrowsize = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*120))) + int(1) edgelabel = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate))) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') else: linewidth = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*20))) + int(1) arrowsize = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate*120))) + int(1) edgelabel = abs(int(((float(openmkm_reaction_netrate[j]))*100)/(normalization_rate))) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') if count > 0: color1 = get_color(openmkm_species_surf_cov[i],openmkm_species_phase[i]) f.write('"' + str(openmkm_species_list[i]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str("Species\\") + str(i) + '.svg",shape=plaintext];\n') def normalization_local(fname): print('now doing local consumption visualization') # print('IN here',openmkm_species_list) for i in range(len(openmkm_species_list)): if openmkm_species_consider[i] == True: count = 0 #Keeping a count of how many reactions with reaction rate more than cutoff are present for j in range(len(openmkm_reaction_strings)): #First checking if a reaction exists above the cutoff rate if float(stoich_matrix[i][j] != 0.0): if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): count = count + 1 #So if the species is being consumed in the reaction if float(rpa_local[i][j]) <= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: # print(openmkm_species_list[i],prod) #check if the product needs to be considered for visualization if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') if float(rpa_local[i][j]) >= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: # print(openmkm_species_list[i], prod, 'new') if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(edgelabel) + '%"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') if count > 0: #If surface coverages provided change fillcolor based on the legend color1 = get_color(openmkm_species_surf_cov[i],openmkm_species_phase[i]) f.write('"' + str(openmkm_species_list[i]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str("Species\\") + str(i) + '.svg",shape=plaintext];\n') openmkm_species_list = [] openmkm_species_phase = [] openmkm_species_surf_cov = [] openmkm_species_composition = [] openmkm_reaction_strings = [] openmkm_reaction_reactants = [] openmkm_reaction_products = [] openmkm_reaction_fwdrate = [] openmkm_reaction_revrate = [] openmkm_reaction_netrate = [] openmkm_reaction_pei = [] #Inputs from OpenMKM or any other chemical kinetics software InitialReactant = 'NH3' Equilibrium_Tolerance = 0.05 Reaction_Rate_Cutoff = 1.0E-09 Elemental_List = ['N', 'H'] Elements_Available = [] Rank_Sep = 0.25 Node_Sep = 0.25 Normalization_requested = 1 Equilibrium_Upper = 0.5 + Equilibrium_Tolerance Equilibrium_Lower = 0.5 - Equilibrium_Tolerance #Implementing pandas to read in the header and initialize the element fields #delimiter = '\t' #pandas puts in the delimited \t between spaces and #so it needs to be removed when reading in the file #Populating all elements provided in input species file input_species = pd.read_csv('species_comp.out', header = 0, delimiter = '\t') for col in input_species.columns: if col != 'SpeciesName' and col != 'Phase': Elements_Available.append(col) #print(Elements_Available) #Populating species datastructures with open('species_comp.out', 'r') as f_ptr: next(f_ptr) for line in f_ptr: line = line.replace('\n', '') words = line.split() openmkm_species_list.append(words[0]) openmkm_species_phase.append(words[1]) openmkm_species_surf_cov.append(words[2]) composition = {} for j in range(3,len(words)): if int(words[j]) > 0: composition[Elements_Available[j-3]] = words[j] openmkm_species_composition.append(composition) #print('Species stored') #Populating reactions datastructures #print(openmkm_species_surf_cov) #reactions_read('rates_ss.out') with open('rates_ss.out') as f_ptr: for i in range(2): f_ptr.readline() for line in f_ptr: line = line.replace('\n', '') words = line.split() openmkm_reaction_fwdrate.append(words[0]) openmkm_reaction_revrate.append(words[1]) openmkm_reaction_netrate.append(words[2]) openmkm_reaction_pei.append(words[3]) openmkm_reaction_strings.append(words[4]) openmkm_species_consider = [] for i in range(len(openmkm_species_list)): # print(openmkm_species_list[i]) consider = False species_dict = openmkm_species_composition[i] # print(species_dict) #Checking if elements in user specified list is present in species for elem in Elemental_List: for elem1 in species_dict.keys(): if elem == elem1: consider = True openmkm_species_consider.append(consider) filename = normalization_type(Normalization_requested) #print(filename) pre, ext = os.path.splitext(filename) #print(pre,ext) filename_svg = str(pre) + str('.svg') #print(filename_svg) print_header(filename, Rank_Sep, Node_Sep) f = open(filename, "a+") #General output specifying number of reactions and species in system TotalReactions = len(openmkm_reaction_netrate) TotalSpecies = len(openmkm_species_list) #print('The networks includes', TotalReactions, ' reactions and ',TotalSpecies, ' species!') #print('Checking if initial reactant is present: ', openmkm_species_list.count(InitialReactant)) MaxReactionRate = max(abs(float(n)) for n in openmkm_reaction_netrate) #print('Maximum reaction rate is: ',MaxReactionRate) stoich_matrix = np.zeros( (TotalSpecies, TotalReactions) ) total_corelations = 0 openmkm_reactions_list = [] with open('reactions.out','r') as f_ptr: for line in f_ptr: line = line.replace('\n', '') reac_prod = line.split('<=>') reactants = reac_prod[0] products = reac_prod[1] reaction_index = openmkm_reaction_strings.index(line) reactant = reactants.split('+') reactant_without_coeff = [] for reac in reactant: reac_split = reac.split('|') if len(reac_split) > 1: reac_stoich = reac_split[0] reac_string = reac_split[1] else: reac_stoich = 1 reac_string = reac_split[0] species_index = openmkm_species_list.index(reac_string) stoich_matrix[species_index][reaction_index] = -int(reac_stoich) reactant_without_coeff.append(reac_string) total_corelations = total_corelations + 1 openmkm_reaction_reactants.append(reactant_without_coeff) product = products.split('+') product_without_coeff = [] for prod in product: prod_split = prod.split('|') if len(prod_split) > 1: prod_stoich = prod_split[0] prod_string = prod_split[1] else: prod_stoich = 1 prod_string = prod_split[0] species_index = openmkm_species_list.index(prod_string) stoich_matrix[species_index][reaction_index] = int(prod_stoich) product_without_coeff.append(prod_string) total_corelations = total_corelations + 1 openmkm_reaction_products.append(product_without_coeff) #print('total number of rows in rpa_local matrix ',total_corelations) reaction_deln = np.sum(stoich_matrix, axis=0) #print('Listing all deln of reactions: ',reaction_deln) initial_reactant_index = openmkm_species_list.index(InitialReactant) #print('Index of initial reactant: ',initial_reactant_index) NetRate = 0.0 #print(stoich_matrix[initial_reactant_index]) for i in range(len(stoich_matrix[initial_reactant_index])): NetRate = NetRate + float(stoich_matrix[initial_reactant_index][i])*float(openmkm_reaction_netrate[i]) #print('NetRate of starting species: ',abs(NetRate),' species: ',InitialReactant) rpa_local = np.zeros( (TotalSpecies, TotalReactions) ) rpa_consumption = [] rpa_production = [] #Calculate the net rate of formation and consumption for each species for i in range(len(openmkm_species_list)): NetConsumptionRate = 0.0 NetProductionRate = 0.0 for j in range(len(openmkm_reaction_strings)): if float(stoich_matrix[i][j]) < 0 and float(openmkm_reaction_netrate[j]) > 0: # => consumption NetConsumptionRate = NetConsumptionRate + float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j]) elif float(stoich_matrix[i][j]) < 0 and float(openmkm_reaction_netrate[j]) < 0: # => production NetProductionRate = NetProductionRate + float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j]) elif float(stoich_matrix[i][j]) > 0 and float(openmkm_reaction_netrate[j]) > 0: # => production NetProductionRate = NetProductionRate + float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j]) elif float(stoich_matrix[i][j]) > 0 and float(openmkm_reaction_netrate[j]) < 0: # => consumption NetConsumptionRate = NetConsumptionRate + float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j]) else: # print('Either species does not participate in reaction or reaction rate is zero') pass rpa_consumption.append(NetConsumptionRate) rpa_production.append(NetProductionRate) #print(rpa_consumption) #print(rpa_production) #Next generate rpa_local for each species for i in range(len(openmkm_species_list)): for j in range(len(openmkm_reaction_strings)): if float(stoich_matrix[i][j]) < 0 and float(openmkm_reaction_netrate[j]) > 0: # => consumption rpa_local[i][j] = float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j])*float(100)/float(abs(rpa_consumption[i])) elif float(stoich_matrix[i][j]) < 0 and float(openmkm_reaction_netrate[j]) < 0: # => production rpa_local[i][j] = float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j])*float(100)/float(abs(rpa_production[i])) elif float(stoich_matrix[i][j]) > 0 and float(openmkm_reaction_netrate[j]) > 0: # => production rpa_local[i][j] = float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j])*float(100)/float(abs(rpa_production[i])) elif float(stoich_matrix[i][j]) > 0 and float(openmkm_reaction_netrate[j]) < 0: # => consumption rpa_local[i][j] = float(openmkm_reaction_netrate[j])*float(stoich_matrix[i][j])*float(100)/float(abs(rpa_consumption[i])) else: # print('Either species does not participate in reaction or reaction rate is zero') pass #print(rpa_local) if Normalization_requested == 1: normalization(filename, MaxReactionRate) elif Normalization_requested == 2: normalization(filename, abs(NetRate)) elif Normalization_requested == 3: normalization_local(filename) else: f.write('This type of normalization does not exist\n') f.write('}') f.close() #Species Visualizations need to be added check_call(['dot','-Tsvg',f.name,'-o',filename_svg]) generate_legend('legend_cov.out') check_call(['dot','-Tsvg','legend_cov.out','-o','legend.svg']) #check_call(['dot','-Tpng','RefinedVisualization.dot','-o','rpa_visualizationMaxRateNormalized1.png']) #check_call(['dot','-Tsvg','RefinedVisualization.dot','-o','rpa_visualizationMaxRateNormalized1.svg']) #dot -Ttiff rpa_visualizationMaxRateNormalized.out -o rpa_visualizationMaxRateNormalized1.svg #(graph,) = pydot.graph_from_dot_file('rpa_visualizationMaxRateNormalized.dot') #graph.write_png('rpa_visualizationMaxRateNormalized1.png') #Check if species folder is present check = os.path.exists("Species") #print(os.path.exists("Species")) if check == False: os.mkdir('Species') #print(os.path.exists("Species")) #Generate species visualizations for i in range(len(openmkm_species_list)): speciescounter = int(i) # print(speciescounter) fname = str('Species/') + str(speciescounter) + str('.txt') fname_svg = str('Species/') + str(speciescounter) + str('.svg') print_header_species(fname, Rank_Sep, Node_Sep) f = open(fname, "a+") color = get_color(openmkm_species_surf_cov[i],openmkm_species_phase[i]) f.write('"' + str(openmkm_species_list[i]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color) + ',URL="' + str(i) + '.svg",shape=plaintext];\n') for j in range(len(openmkm_reaction_strings)): if float(rpa_local[i][j]) <= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): color1 = get_color(openmkm_species_surf_cov[int(openmkm_species_list.index(prod))],openmkm_species_phase[int(openmkm_species_list.index(prod))]) f.write('"' + str(openmkm_species_list[int(openmkm_species_list.index(prod))]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str(int(openmkm_species_list.index(prod))) + '.svg",shape=plaintext];\n') if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(prod) + '"\n') if float(rpa_local[i][j]) <= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) > 0.0: for reac in openmkm_reaction_reactants[j]: if openmkm_species_consider[int(openmkm_species_list.index(reac))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): color1 = get_color(openmkm_species_surf_cov[int(openmkm_species_list.index(reac))],openmkm_species_phase[int(openmkm_species_list.index(reac))]) f.write('"' + str(openmkm_species_list[int(openmkm_species_list.index(reac))]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str(int(openmkm_species_list.index(reac))) + '.svg",shape=plaintext];\n') if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(reac) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) < 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(openmkm_species_list[i]) + '"->"' + str(reac) + '"\n') if float(rpa_local[i][j]) >= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) < 0.0: for prod in openmkm_reaction_products[j]: if openmkm_species_consider[int(openmkm_species_list.index(prod))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): color1 = get_color(openmkm_species_surf_cov[int(openmkm_species_list.index(prod))],openmkm_species_phase[int(openmkm_species_list.index(prod))]) f.write('"' + str(openmkm_species_list[int(openmkm_species_list.index(prod))]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str(int(openmkm_species_list.index(prod))) + '.svg",shape=plaintext];\n') if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(prod) + '"->"' + str(openmkm_species_list[i]) + '"\n') if float(rpa_local[i][j]) >= 0.0 and float(stoich_matrix[i][j]) != 0.0 and float(stoich_matrix[i][j]) > 0.0: for reac in openmkm_reaction_reactants[j]: if openmkm_species_consider[int(openmkm_species_list.index(reac))] == True: if abs(float(openmkm_reaction_netrate[j])) >= float(Reaction_Rate_Cutoff): color1 = get_color(openmkm_species_surf_cov[int(openmkm_species_list.index(reac))],openmkm_species_phase[int(openmkm_species_list.index(reac))]) f.write('"' + str(openmkm_species_list[int(openmkm_species_list.index(reac))]) + '"[shape=rectangle,style=filled,fontsize=35,width=0,height=0,fillcolor=' + str(color1) + ',URL="' + str(int(openmkm_species_list.index(reac))) + '.svg",shape=plaintext];\n') if float(openmkm_reaction_pei[j]) >= Equilibrium_Lower and float(openmkm_reaction_pei[j]) <= Equilibrium_Upper: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=green,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(reac) + '"->"' + str(openmkm_species_list[i]) + '"\n') else: linewidth = abs(int(rpa_local[i][j]/20)) + int(1) arrowsize = abs(int(rpa_local[i][j]/60)) + int(1) edgelabel = abs(int(rpa_local[i][j])) if float(rpa_local[i][j]) > 0.0: f.write('edge[dir="forward",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') else: f.write('edge[dir="none",style="setlinewidth(' + str(linewidth) + ')",color=black,weight=2,arrowsize=' + str(arrowsize) + ',label=" ' + str(j) + ' ' + str(edgelabel) + '% ' + str(abs(float(openmkm_reaction_netrate[j]))) + ' mol/s ' + str(openmkm_reaction_pei[j]) + '"];\n') f.write('"' + str(reac) + '"->"' + str(openmkm_species_list[i]) + '"\n') f.write('}') f.close() check_call(['dot','-Tsvg',f.name,'-o',fname_svg])