branching and hyperbranching

.gitignore

@@ -139,4 +139,5 @@ tests/output/
 
 polyconf_examples/PEI_linear_*
 polyconf_examples/PEI_branched_*
+polyconf_examples/PEI_randombranched_*
 polyconf_examples/PMMA_*

polyconf_examples/04_build_randomly_branched_PEI.py

@@ -0,0 +1,174 @@
+#!/usr/bin/env python
+
+from polyconf.Monomer import Monomer
+from polyconf.Polymer import Polymer
+from polyconf.PDB import PDB
+
+import random
+
+random.seed(10) # by setting the random seed we can ensure the random conformations are replicable.  
+
+###############################################################################
+#                                                                             #
+#  Branched PEI                                                               #
+#                                                                             #
+###############################################################################
+
+# In this example, we will follow the process used to generate a linear PEI, 128 monomers long
+
+# but then we will add branches every 5 monomers 
+
+PEI_branched=Polymer(Monomer('PEI_start.pdb')) 
+
+adds=126 
+
+for i in range (0,adds):
+    PEI_branched.extend( 
+            Monomer('PEI_monomer.pdb'),
+            n=PEI_branched.maxresid(), 
+            nn=PEI_branched.newresid(),
+            names=dict(P='C1',Q='CX',R='NX',S='N1'), 
+            joins=[('N1','C1')],
+            ) 
+
+# add the final monomer
+
+PEI_branched.extend(
+            Monomer('PEI_end.pdb'),
+            n=PEI_branched.maxresid(), 
+            nn=PEI_branched.newresid(),
+            names=dict(P='C1',Q='CX',R='NX',S='N1'),
+            joins=[('N1','C1')],
+            )
+
+
+# now we have a linear chain, we can add our branches
+# for this example each branch will be five units long
+
+# we want to attach our branches randomly
+# this means we want to choose randomly from a list of possible attachment points 
+# and that means we need a way to get a list of possible attachment points
+
+
+
+def get_hardpoints(polymer,selstr):
+        # this function will take a polymer and a selection string, and return a list of all resids that match the sleection string
+        # we will ask for a list of all monomers where the residue name is ELNR, which matches the 126 middle monomers
+        # where the monomer contains an H1 atom
+        # this will exclude monomers where the H1 atom has been removed because we have added a branch
+        atomlist = polymer.select_atoms(selstr)
+        return(list(atomlist.atoms.resids))
+
+hardpoints = get_hardpoints(PEI_branched,"resname ELNR and name H1")
+
+
+
+branch_bonds=[] # make a list to keep track of all of our branches
+
+
+
+branch=1 # counter to label the branches
+
+while branch <= 10:
+
+    hardpoints = get_hardpoints(PEI_branched,"resname ELNR and name H1") # generate a list of all possible attachment sites
+
+    n = random.sample(hardpoints,1)[0] # chose one of the available hardpoints to branch out of
+    # if we have alreaday created any branches, this list will include the four ELNR monomers from those branches
+
+
+    branch_resid=PEI_branched.newresid()
+    branch_bonds+=[(n,branch_resid)]
+    # add the first monomer of the branch
+    PEI_branched.extend( 
+            Monomer('PEI_monomer.pdb'),
+            n=n, 
+            nn=branch_resid,
+            names=dict(P='C1',Q='H1',R='NX',S='N1'),
+            joins=[('N1','C1')],
+            beta=branch # we are using beta factors to label our branches
+            ) 
+
+    PEI_branched.renamer(n,'H1')
+
+    # add the next three monomers of the branch
+
+    for _ in range(0,3):
+            PEI_branched.extend(
+                Monomer('PEI_monomer.pdb'),
+                n=PEI_branched.maxresid(), 
+                nn=PEI_branched.newresid(),
+                names=dict(P='C1',Q='CX',R='NX',S='N1'),
+                joins=[('N1','C1')],
+                beta=branch
+
+                )
+
+    # add the fifth monomer of the branch
+
+    PEI_branched.extend(
+            Monomer('PEI_end.pdb'),
+            n=PEI_branched.maxresid(), 
+            nn=PEI_branched.newresid(),
+            names=dict(P='C1',Q='CX',R='NX',S='N1'),
+            joins=[('N1','C1')],
+            beta=branch
+
+            )
+
+    branch+=1 # change the branch number for the next branch
+
+# Generate lists of dihedrals 
+
+
+
+mainchain_NC= PEI_branched.gen_pairlist(J='N1',K='C1',first_resid=1,last_resid=127,mult=3,same_res=False) 
+mainchain_alkanes= PEI_branched.gen_pairlist(J='C1',K='C2',first_resid=1,last_resid=128,mult=3)
+
+# we are going to manually generate a dihedral list for the connections between the main chain and the branches
+branch_dihedrals=[{'J': 'N1', 'J_resid': i, 'K': 'C1', 'K_resid': j, 'mult': 3} for i,j in branch_bonds]
+
+branch_alkanes=PEI_branched.gen_pairlist(J='C1',K='C2',first_resid=129,last_resid=PEI_branched.maxresid(),mult=3)
+
+for conf in range(1,6): 
+
+    print(f'generating conformation {conf}')
+    newconf=PEI_branched.copy() # make a duplicate of the original polymer
+
+    print(f'shuffling backbone')
+
+    newconf.shuffler(mainchain_NC)  # as before, randomise the NC dihedrals to generate a random conformation 
+
+    print(f'solving backbone')
+
+    newconf.dihedral_solver(mainchain_alkanes,dummies='CX NX X*',cutoff=0.9) # solve by rotating the alkane dihedrals
+
+    print(f'shuffling branchpoints')
+
+    newconf.shuffler(branch_dihedrals) 
+
+    print(f'solving branchpoints')
+
+    newconf.dihedral_solver(branch_dihedrals,dummies='CX NX X*',cutoff=0.9,backwards_only=False) 
+
+    print(f'solving branches')
+
+    newconf.dihedral_solver(branch_alkanes,dummies='CX NX X*',cutoff=0.9,backwards_only=False) 
+
+
+
+    Saver = PDB(newconf)
+    Saver.cleanup()
+    Saver.save(dummies='CX NX X1',fname=f'PEI_randombranched_conformation_{conf}') 
+
+# This should produce an ensemble of five starting conformations
+
+# when you visualize these, try colouring by beta factor again
+
+# in particular, look at branches where the beta factor is 1, 4 or 6
+
+# you can see that branch 1 extends from the original PEI backbone
+# but then, branch 4 extends off of branch 1
+# and then, branch 6 extends off of branch 4
+
+# by designing a branching algorithm suitable for your polymer of interest, you can use this type of approach to generate recursively branched or hyperbranched polymers

polyconf_examples/README.md

@@ -8,9 +8,7 @@ These tutorials are not intended as a primer on polymer design or stereochemistr
 
 ## Tutorial 01:  Building a linear polyethyleneimine polymer conformation using PolyConf
 
-**TODO:** extend based on rewrite
-
-The first tutorial is an example of how to make a linear 128 unit polyethyleneimine chain using *PolyConf*. It is contained in the file `01_build_PEI_linear.py`
+The first tutorial is an example of how to make a linear 128 unit polyethyleneimine chain using *PolyConf*. 
 
 This tutorial showcases several approaches, including:
 
@@ -19,7 +17,27 @@ This tutorial showcases several approaches, including:
 * using `dihedral_solver()` to adjust the polyethyleneimine conformation
 * using `shuffler()` and `dihedral_solver()` to generate a random polyethyleneimine conformation
 
-This tutorial is all contained in a single script.  There are detailed comments discussing the strategies used, and the limitations of the resulting conformations.  You are encouraged to examine the conformations generated by this script in [VMD](https://www.ks.uiuc.edu/Research/vmd/), or another visualization tool of your choice, to understand the strengths and limitations of the different strategies.
+This tutorial is contained in a series of scripts.  There are detailed comments discussing the strategies used, and the limitations of the resulting conformations.  You are encouraged to examine the conformations generated by this script in [VMD](https://www.ks.uiuc.edu/Research/vmd/), or another visualization tool of your choice, to understand the strengths and limitations of the different strategies.
+
+### Building a simple linear polymer
+
+`01a_build_PEI_with_extend.py` shows an example of building a polyethyleneimine 128-mer using extend, sand saving the resulting conformation.
+
+### Building a simple linear polymer fit to a linear vector `extend(linearise=True)`
+
+`01b_build_PEI_with_linear_extend.py` shows an example of building a polyethyleneimine 128-mer where the monomers are extended along a specified vector, and showcases the strengths and weaknesses of this approach
+
+### Solving a polymer conformation
+
+`01c_build_PEI_then_solve.py` shows how to use `dihedral_solver()` to convert the conformation from tutorial 01a into one that is more physically reasonable.
+
+### Generating a random polymer conformation
+
+`01d_build_PEI_conformation.py` shows how to use `shuffler()` and `dihedral_solver()` to convert the highly ordered conformation from tutorial 01a into a random coil that is still physically reasonable.
+
+### Generating an ensemble of random polymer conformation
+
+`01e_build_PEI_conformation_ensemble.py` shows how to use the methods from tutorial 01d to generate an ensemble of physically reasonable polymer conformations for replicate molecular dynamics simulations.
 
 ## Tutorial 02:  Strategies for building conformations of polymers with different types of tacticity, or copolymers with different monomer distributions
 
@@ -57,10 +75,22 @@ This example is one solution to resolve the difficulties showcased in `02c_build
 
 Again, this type of approach could also be used to generate a [statistical copolymer](https://en.wikipedia.org/wiki/Copolymer#Statistical_copolymers).
 
-## Tutorial 03:  Building a branched polyethyleneimine polymer conformation using PolyConf
+## Tutorial 03:  Building a regularly branched polyethyleneimine polymer conformation using PolyConf
+
+This tutorial showcases how to generate a starting conformation for a branched polyethyleneimine chain with regularly spaced branches. It is contained in the script `03_build_branched_PEI.py`
+
+First, the script generates a linear polyethyleneimine 128-mer
+Then, the linear 128-mer is extended by adding 3mer branches on every 5th monomer
+
+Finally, a set of random conformations are generated.  This process includes a demostration of how to manually generate a list of dihedrals.
+
+## Tutorial 04:  Building a hyperbranched polyethyleneimine polymer conformation with random branchpoints
+
+This tutorial showcases how to generate a starting conformation for a branched polyethyleneimine chain with randomly spaced branches. It is contained in the script `04_build_randomly_branched_PEI.py`
 
-**TODO:** Ada to write description
+First, the script generates a linear polyethyleneimine 128-mer
+Then, the linear 128-mer is extended by adding 5mer branches on ten randomly selected monomers.
 
-## Tutorial 04:  Building a polyethyleneimine dendrimer conformation using PolyConf
+The method used demonstrates random branch sites, and hyperbranching.  Branches can extend off of the original 128mer, but they can also extend of of other branches.
 
-To be written
+Finally, a set of random conformations are generated.