Files
nmrproject/mod/butadien/butadienCycleExclusion.py
T

277 lines
8.0 KiB
Python

config.ilp.solver="CPLEX"
#import mod
#from mod import *
import networkx as nx
# import Graph class from graph.py
from graph import GraphObj
#ethylen = Graph.fromGMLString(
"""graph [
node [ id 0 label "H" ]
node [ id 1 label "H" ]
node [ id 2 label "C" ]
node [ id 3 label "C" ]
node [ id 4 label "H" ]
node [ id 5 label "H" ]
edge [ source 0 target 2 label "-" ]
edge [ source 1 target 2 label "-" ]
edge [ source 2 target 3 label "=" ]
edge [ source 3 target 4 label "-" ]
edge [ source 3 target 5 label "-" ]
]"""
#, name="Ethylen")
butadien = Graph.fromGMLString(
"""graph [
node [ id 0 label "C" ]
node [ id 1 label "C" ]
node [ id 2 label "C" ]
node [ id 3 label "C" ]
node [ id 4 label "H" ]
node [ id 5 label "H" ]
node [ id 6 label "H" ]
node [ id 7 label "H" ]
node [ id 8 label "H" ]
node [ id 9 label "H" ]
edge [ source 0 target 1 label "=" ]
edge [ source 1 target 2 label "-" ]
edge [ source 2 target 3 label "=" ]
edge [ source 0 target 4 label "-" ]
edge [ source 0 target 5 label "-" ]
edge [ source 1 target 6 label "-" ]
edge [ source 2 target 7 label "-" ]
edge [ source 3 target 8 label "-" ]
edge [ source 3 target 9 label "-" ]
]"""
, name="Butadien")
#benzene = Graph.fromSMILES('c1=cc=cc=c1')
#benzene = Graph.fromGMLString(
"""graph [
node [ id 0 label "C" ]
node [ id 1 label "H" ]
node [ id 2 label "C" ]
node [ id 3 label "H" ]
node [ id 4 label "C" ]
node [ id 5 label "H" ]
node [ id 6 label "C" ]
node [ id 7 label "H" ]
node [ id 8 label "C" ]
node [ id 9 label "H" ]
node [ id 10 label "C" ]
node [ id 11 label "H" ]
edge [ source 0 target 1 label "-" ]
edge [ source 0 target 2 label "=" ]
edge [ source 2 target 3 label "-" ]
edge [ source 2 target 4 label "-" ]
edge [ source 4 target 5 label "-" ]
edge [ source 4 target 6 label "=" ]
edge [ source 6 target 7 label "-" ]
edge [ source 6 target 8 label "-" ]
edge [ source 8 target 9 label "-" ]
edge [ source 8 target 10 label "=" ]
edge [ source 10 target 11 label "-" ]
edge [ source 10 target 0 label "-" ]
]"""
#, name="Benzene")
#pentadien = Graph.fromGMLString(
"""graph
[
node [ id 0 label "C" ]
node [ id 1 label "C" ]
node [ id 2 label "C" ]
node [ id 3 label "C" ]
node [ id 4 label "H" ]
node [ id 5 label "H" ]
node [ id 6 label "H" ]
node [ id 7 label "H" ]
node [ id 8 label "H" ]
node [ id 9 label "C" ]
node [ id 10 label "H" ]
node [ id 11 label "H" ]
node [ id 12 label "H" ]
edge [ source 0 target 1 label "=" ]
edge [ source 1 target 2 label "-" ]
edge [ source 2 target 3 label "=" ]
edge [ source 0 target 4 label "-" ]
edge [ source 0 target 5 label "-" ]
edge [ source 1 target 6 label "-" ]
edge [ source 2 target 7 label "-" ]
edge [ source 3 target 8 label "-" ]
edge [ source 3 target 9 label "-" ]
edge [ source 9 target 10 label "-" ]
edge [ source 9 target 11 label "-" ]
edge [ source 9 target 12 label "-" ]
]
"""
#, name="Pentadien")
restswap = Rule.fromGMLString(
"""rule [
left [
edge [ source 1 target 2 label "=" ]
edge [ source 3 target 4 label "=" ]
]
context [
node [ id 1 label "C" ]
node [ id 2 label "C"]
node [ id 3 label "C"]
node [ id 4 label "C"]
]
right [
edge [ source 1 target 3 label "=" ]
edge [ source 2 target 4 label "=" ]
]
]"""
)
dielsalder = Rule.fromGMLString(
"""rule [
left [
edge [ source 1 target 2 label "=" ]
edge [ source 2 target 3 label "-" ]
edge [ source 3 target 4 label "=" ]
edge [ source 5 target 6 label "=" ]
]
context [
node [ id 1 label "C" ]
node [ id 2 label "C"]
node [ id 3 label "C"]
node [ id 4 label "C"]
node [ id 5 label "C"]
node [ id 6 label "C"]
]
right [
edge [ source 1 target 2 label "-" ]
edge [ source 2 target 3 label "=" ]
edge [ source 3 target 4 label "-" ]
edge [ source 4 target 5 label "-" ]
edge [ source 5 target 6 label "-" ]
edge [ source 6 target 1 label "-" ]
]
]"""
)
def cyclesizes(g):
nxGraph = GraphObj(g).nx_graph
#Restriction to E/Z? Benzen can't be created or consumed by the reaction due to steric difficulties
#Chain restrictions: (No steric reason)
for nodestart in nxGraph.nodes:
for nodeend in nxGraph.nodes:
if nodeend != nodestart:
paths = sorted(list(nx.all_shortest_paths(nxGraph, nodestart, nodeend)))
for path in paths:
if len(path) > 10: #First node included in length
return True
#Cycle restrictions:
cycles = sorted(list(nx.chordless_cycles(nxGraph)))
cycle_lengths = [len(x) for x in cycles]
if cycles == []:
return False
#One Atom can't be in four different cycles
for vertice in nxGraph.nodes:
counter = 0
for cycle in cycles:
if vertice in cycle:
counter += 1
if counter >= 4:
return True
#One Cycle can't overlapp with another on over 2 Connection points
for cycle in cycles:
for cycleref in cycles:
if cycle != cycleref and len(list(set(cycle) & set(cycleref))) >= 3:
return True
#Only cordless cycles of length 5,6 and 7 are acceptable
if min(cycle_lengths) >= 6 and max(cycle_lengths) <=6:
return False
return True
#Disallows Allenes (Two Doublebonds on same carbon)
def doubledoublebond(g):
for vertice in g.vertices:
if (vertice.stringLabel == "C" and vertice.degree == 2): # and vertice.edges.label == ["=", "="]
for edge in vertice.incidentEdges:
print(type(edge.bondType))
bonds = [type(edge.bondType) for edge in vertice.incidentEdges]
if (bonds[0] == bonds[1] and len(bonds) == 2):
return True
return False
def restriction(dg):
for a in dg.right:
if a.vLabelCount("C") > 10:
return False
if doubledoublebond(a):
return False
if cyclesizes(a):
return False
return True
flowPrinter = FlowPrinter()
flowPrinter.printUnfiltered = False
postSection("Loaded Graphs")
for a in inputGraphs:
a.print()
postSection("Loaded Rules")
for a in inputRules:
a.print()
ethylen = Graph.fromSMILES('C=C')
goal = Graph.fromSMILES('C1C=CC=CC=1')
dg = DG(graphDatabase=inputGraphs)
dg.build().execute(
addSubset(inputGraphs)
>> rightPredicate[
restriction
](
repeat(revive(inputRules)) #Revive not necessary
)
)
print([str(a.graph.name) for a in dg.vertices])
print([str(a.graph.smiles) for a in dg.vertices])
dg.print()
postSection("Product Graphs")
for a in dg.vertices:
a.graph.print()
flow = Flow(dg)
flow.addSource(butadien)
flow.addSink(ethylen)
flow.addSink(goal)
flow.addConstraint(inFlow[butadien] >= 2)
#flow.addConstraint(outFlow[ethylen] == 1)
flow.addConstraint(outFlow[goal] == 1)
flow.findSolutions(maxNumSolutions=3)
flow.solutions.list()
flow.solutions.print(flowPrinter)
sys.exit(0)
rc = rcEvaluator(inputRules)
for dRef in dg.derivations:
der = dRef.derivation
educt = rcId(der.left[0])
for i in range(1, len(der.left)):
educt = educt *rcParallel* rcId(der.left[i])
product = rcId(der.right[0])
for i in range(1, len(der.right)):
product = product *rcParallel* rcId(der.right[i])
rcExp = educt *rcSuper(allowPartial=False)* der.rule *rcSuper(allowPartial=False)* product
res = rc.eval(rcExp)
dRef.print()
for a in res:
a.print()
a.printGML()