diff --git a/HTPolyNet/molecule.py b/HTPolyNet/molecule.py index 8ababdb..e031f6e 100644 --- a/HTPolyNet/molecule.py +++ b/HTPolyNet/molecule.py @@ -127,7 +127,7 @@ def __init__(self,name='',generator:Reaction=None,origin:str=None): self.bond_templates:BondTemplateList=[] self.symmetry_relateds=[] self.stereocenters=[] # list of atomnames - self.stereoisomers:dict(str,Molecule)={} + self.stereoisomers={} self.nconformers=0 self.conformers_dict={} self.conformers=[] # just a list of gro file basenames @@ -890,7 +890,7 @@ def transrot(self,at_idx,at_resid,from_idx,from_resid,connected_resids=[]): overall_maximum=(-1.e9,-1,-1) coord_trials={} for myH,myHnm in myHpartners.items(): # keys are globalIdx's, values are names - coord_trials[myH]:dict[TopoCoord]={} + coord_trials[myH]={} Rh=TC.get_R(myH) logger.debug(f' Rh {myH} {Rh} {Rh.dtype}') Rih=Ri-Rh diff --git a/HTPolyNet/software.py b/HTPolyNet/software.py index 8178935..178d3fe 100644 --- a/HTPolyNet/software.py +++ b/HTPolyNet/software.py @@ -10,7 +10,6 @@ import logging import os from HTPolyNet.stringthings import my_logger -# from GPUtil import getGPUs logger=logging.getLogger(__name__) class Software: diff --git a/HTPolyNet/topocoord.py b/HTPolyNet/topocoord.py index bde0581..1436241 100644 --- a/HTPolyNet/topocoord.py +++ b/HTPolyNet/topocoord.py @@ -21,6 +21,7 @@ from HTPolyNet.matrix4 import Matrix4 from HTPolyNet.gromacs import grompp_and_mdrun,mdp_get, mdp_modify, gmx_energy_trace import HTPolyNet.projectfilesystem as pfs +# import HTPolyNet.polybonds as pb logger=logging.getLogger(__name__) @@ -31,9 +32,9 @@ class BTRC(Enum): :type Enum: class """ passed = 0 - failed_pierced_ring = 1 # does this bond-candidate pierce a ring? - failed_shortcircuit = 2 # does this bond-candidate result in short-circuit? - failed_bondcycle = 3 # does this bond-candidate create a bondcycle on its own? + failed_pierced_ring = 1 # candidate bond pierces a ring + failed_shortcircuit = 2 # candidate bond creates a short-circuit + # failed_bondcycle = 3 # Bondcycles are handled collectively unset = 99 class TopoCoord: @@ -65,6 +66,7 @@ def __init__(self,topfilename='',tpxfilename='',grofilename='',grxfilename='',mo self.grxattr=[] self.idx_lists={} self.idx_lists['bondchain']=[] + # self.polybondchainmanager=pb.PolyBondChainManager() if grofilename!='': self.read_gro(grofilename,wrap_coords=wrap_coords) else: @@ -1241,8 +1243,8 @@ def bondtest(self,b,pbc=[1,1,1],show_piercings=True): # j=int(j) if self.makes_shortcircuit(i,j): return BTRC.failed_shortcircuit,0 - if self.makes_bondcycle(i,j): - return BTRC.failed_bondcycle,0 + # if self.makes_bondcycle(i,j): + # return BTRC.failed_bondcycle,0 if self.pierces_ring(i,j,pbc=pbc,show_piercings=show_piercings): return BTRC.failed_pierced_ring,0 # logger.debug(f'passed: {i:>7d} {j:>7d} {rij:>6.3f} nm') @@ -1461,84 +1463,148 @@ def bondchainlist_update(self,new_bond_recs,msg=''): cnms.append([self.get_gro_attribute_by_attributes('atomName',{'globalIdx':x}) for x in c]) # logger.debug(f'post {msg} chains {self.idx_lists["bondchain"]} {cnms}') - def makes_bondcycle(self,aidx,bidx): - """makes_bondcycle checks the current 'bondchain' index lists to see if a bond between aidx and bidx (global atom indices) would generate a C-C' bondcycle, which they do if they have the same 'bondchain' attribute, shorter than the minimum allowable bondcycle length, if this is not -1 + # def makes_bondcycle(self,aidx,bidx): + # """makes_bondcycle checks the current 'bondchain' index lists to see if a bond between aidx and bidx (global atom indices) would generate a C-C' bondcycle, which they do if they have the same 'bondchain' attribute, shorter than the minimum allowable bondcycle length, if this is not 0 - :param aidx: global index of an atom - :type aidx: int - :param bidx: global index of another atom - :type bidx: int - :return: True if a bond betwen aidx and bidx would form a C-C' bondcycle - :rtype: bool - """ - # is there a bondchain with aidx as head and bidx as tail, or vice versa? - for c in self.idx_lists['bondchain']: - if (aidx==c[0] and bidx==c[-1]) or (aidx==c[-1] and bidx==c[0]): - if self.min_bondcycle_length==-1 or len(c) Cycle length is {bl.len()} reacted C=C bonds') + if self.min_bondcycle_length<=0 or bl.len() limit is {self.min_bondcycle_length}, so we will break the longest added bond') + logger.debug(f'-> bond index {longest_bond_index} is the longest added bond in the cycle') + new_bdf.loc[longest_bond_index,'remove-to-uncyclize']=True + + # if rdict.get('makescycle',False): + # chains_of_bonds[i]=rdict["cyclizedchain"] + # cyclized_chains.append(rdict["cyclizedchain"]) + # else: + # if rdict.get('retainedchain',-1)>-1: + # # a chain merger happened + # chains_of_bonds[i]=rdict["retainedchain"] + # for k in chains_of_bonds.keys(): + # if chains_of_bonds[k]==rdict["removedchain"]: + # chains_of_bonds[k]=rdict["retainedchain"] + # else: + # # no chain merger happened; do nothing + # pass + # # for each bondcycle, identify the new bonds that collectively cause the cyclization + # cyclizing_bonds={} + # for k,v in chains_of_bonds.items(): + # if v in cyclized_chains: + # if not v in cyclizing_bonds: + # cyclizing_bonds[v]=[] + # cyclizing_bonds[v].append(k) + + # for cycleidx,newbonds in cyclizing_bonds.items(): + # cyclelength=len(chains[cycleidx]) + # # # logger.debug(f'bond {i} {r.ai}-{r.aj} ({r.reactantName}) adds to bondchain {chain_idx}') + # # chains_of_bonds.append(chain_idx) + # # if cyclized and not chain_idx in cyclized_chains: + # # cyclized_chains[chain_idx]=[] + # # assert len(chains_of_bonds)==bdf.shape[0] + # for i,r in bdf.iterrows(): + # # logger.debug(f'compiling bonds of cyclized chains: bond {i}') + # cidx=chains_of_bonds[i] + # # logger.debug(f'compiling bonds of cyclized chains: cidx {chains_of_bonds[i]}') + # if cidx in cyclized_chains: + # cyclized_chains[cidx].append(i) + # for k,v in cyclized_chains.items(): + # if self.min_bondcycle_length<=0 or len(v)