netcdf traj update

Liu-group · Feb 22, 2022 · d2bb316 · d2bb316
1 parent 65cf228
commit d2bb316
Show file tree

Hide file tree

Showing 8 changed files with 349 additions and 14 deletions.
diff --git a/autosolvate/__main__.py b/autosolvate/__main__.py
@@ -42,7 +42,7 @@ def main(args=None):
         startmd(args[1:])
     elif args[0] == 'clustergen':
         print('AutoSolvate is starting in command line mode!')
-        print('Running the module to extract solvated cluster (sphere) from MD trajectories of solvent box.')
+        print('Running the module to extract microsolvated clusters from MD trajectories with solvent box.')
         startclustergen(args[1:])
     else:
         print('Invalid syntax for AutoSolvate command line interface.')

diff --git a/autosolvate/generatetrajs.py b/autosolvate/generatetrajs.py
@@ -287,8 +287,7 @@ def writeQMMMMinInput(stepsqmmmmin=250):
         f.write("gpr QMMM Min\n")
         f.write(" &cntrl\n")
         f.write("  imin   = 1,\n")
-        f.write("  maxcyc = 2000,\n")
-        f.write("  ncyc   = 1000,\n")
+        f.write("  maxcyc = "+str(stepsqmmmmin)+",\n")
         f.write("  irest  = 0, ! 0- new simulation 1- restart\n")
         f.write("  ntx    = 1, ! 1-read in coordinates, but not velocity, 5-both\n")
         f.write("  cut    = 8.0,\n")
@@ -342,8 +341,11 @@ def writeQMMMInput(temperature=300, charge=0, stepsqmmm=250, infilename='qmmmhea
         f.write("gpr QMMM "+infilename+"\n")
         f.write(" &cntrl\n")
         f.write("  imin   = 0,\n")
-        f.write("  irest  = 1, ! 0- new simulation 1- restart\n")
-        f.write("  ntx    = 5, ! 1-read in coordinates, but not velocity, 5-both\n")
+        if infilename=='qmmmheat.in':
+          f.write("  ntx    = 1, ! 1-read in coordinates, but not velocity, 5-both\n")
+        else:
+          f.write("  ntx    = 5, ! 1-read in coordinates, but not velocity, 5-both\n")
+          f.write("  irest  = 1, ! 0- new simulation 1- restart\n")
         f.write("  cut    = 8.0,\n")
         f.write("  ig     = -1, !random seed\n")
         f.write("  ntc    = 2, ntf    = 2, !Shake is used for solvent\n")

diff --git a/docs/_data/naphthalene_radical.xyz b/docs/_data/naphthalene_radical.xyz
@@ -0,0 +1,20 @@
+18
+naphthalene radical
+C     2.4584929186   -0.6980401434    0.0000208854
+C     2.4584830542    0.6980208281    0.0000273558
+C     1.2392834454    1.4064616303    0.0000201346
+C    -0.0000127820    0.7187236077    0.0000072068
+C    -1.2393189424    1.4064428097    0.0000009527
+C    -2.4585398474    0.6980627613   -0.0000119130
+C    -2.4584830245   -0.6980052848   -0.0000182990
+C    -1.2392726206   -1.4064393494   -0.0000121035
+C     0.0000166810   -0.7186826023    0.0000005670
+C     1.2392855074   -1.4064696461    0.0000073561
+H     1.2470800358   -2.4919577916    0.0000018617
+H    -1.2470920207   -2.4919275393   -0.0000168975
+H    -3.3951566422   -1.2429180456   -0.0000277271
+H    -3.3952681112    1.2428765068   -0.0000168560
+H    -1.2469606339    2.4919363439    0.0000063915
+H     1.2471333000    2.4919523490    0.0000239494
+H     3.3951743890    1.2429028846    0.0000376679
+H     3.3951863936   -1.2429173191    0.0000261673
diff --git a/docs/_data/water_solvated-cutoutn-0.xyz b/docs/_data/water_solvated-cutoutn-0.xyz
diff --git a/docs/_data/water_solvated-qmmmheat.netcdf b/docs/_data/water_solvated-qmmmheat.netcdf
diff --git a/docs/_data/water_solvated-qmmmnvt.netcdf b/docs/_data/water_solvated-qmmmnvt.netcdf
diff --git a/docs/_data/water_solvated.prmtop b/docs/_data/water_solvated.prmtop
@@ -1,11 +1,11 @@
-%VERSION  VERSION_STAMP = V0001.000  DATE = 01/01/22  15:59:46                  
+%VERSION  VERSION_STAMP = V0001.000  DATE = 02/19/22  11:39:06                  
 %FLAG TITLE                                                                     
 %FORMAT(20a4)                                                                   
 SLU                                                                             
 %FLAG POINTERS                                                                  
 %FORMAT(10I8)                                                                   
    24003       4   23993      11      16      14      34      20       0       0
-   32069    7996      11      14      20       4       2       2       4       1
+   32069    7996      11      14      20       4       2       2       4       0
        0       0       0       0       0       0       0       1      18       0
        0
 %FLAG ATOM_NAME                                                                 
@@ -18029,7 +18029,7 @@ O   H1  H2
 %FLAG NONBONDED_PARM_INDEX                                                      
 %FORMAT(10I8)                                                                   
        1       2       4       7       2       3       5       8       4       5
-       6      -1       7       8      -1      10
+       6       9       7       8       9      10
 %FLAG RESIDUE_LABEL                                                             
 %FORMAT(20a4)                                                                   
 SLU WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT 
@@ -29737,13 +29737,13 @@ WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT WAT
        0   23999   24000   24000       0   24002   24003   24003       0
 %FLAG HBOND_ACOEF                                                               
 %FORMAT(5E16.8)                                                                 
-  0.00000000E+00
+
 %FLAG HBOND_BCOEF                                                               
 %FORMAT(5E16.8)                                                                 
-  0.00000000E+00
+
 %FLAG HBCUT                                                                     
 %FORMAT(5E16.8)                                                                 
-  0.00000000E+00
+
 %FLAG AMBER_ATOM_TYPE                                                           
 %FORMAT(20a4)                                                                   
 ca  ca  ca  ca  ca  ca  ca  ca  ca  ca  ha  ha  ha  ha  ha  ha  ha  ha  OW  HW  

diff --git a/docs/tutorial.rst b/docs/tutorial.rst
@@ -57,6 +57,16 @@ Once you have AutoSolvate and all dependencies installed you will need the solut
        H     3.3951743890    1.2429028846    0.0000376679
        H     3.3951863936   -1.2429173191    0.0000261673
 
+
+.. note::
+
+  You can download both xyz files here:
+
+  :download:`naphthalene_neutral.xyz <_data/naphthalene_neutral.xyz>`
+
+  :download:`naphthalene_radical.xyz <_data/naphthalene_radical.xyz>`
+
+
 Now that you have the structures, make a directory for each example. We will start with the neutral molecule. 
 
 Example 1: Naphthalene in Water
@@ -309,7 +319,7 @@ To extract the cluster from the final QMMM results, use the following command:
 
   If you were not able to run the QMMM simulation above, you can download the QM/MM trajectory here:
 
-  :download:`water_solvated-qmmmnvt <_data/water_solvated-qmmmnvt.netcdf>`
+  :download:`water_solvated-qmmmnvt.netcdf <_data/water_solvated-qmmmnvt.netcdf>`
 
 The .prmtop and .netcdf filenames are required, but Autosolvate will use the default values of 0 for the starting frame, 100 for the extraction interval, and a cutout size of 4 Å.
 
@@ -333,14 +343,21 @@ The output of this command will be the cartesian coordinates of the microsolvate
 .. image:: _images/tutorial4_3.jpg
    :width: 400
 
-Running the above command only generates one xyz file because we only did 10 steps of the QMMM NVT in our example mdrun, and we asked for a cluster from every ten frames. However, if we extract every step (with option `-i 1`), then we will get 10 coordinate files. 
+Running the above command only generates one xyz file because we only did 100 steps of the QMMM NVT in our example mdrun, and we asked for a cluster from every hundred frames. However, if we extract every 10 steps (with option `-i 10`), then we will get 10 coordinate files. We can increase the solvent shell size to 6 Å with `-s 6`. 
 
-``autosolvate clustergen -f water_solvated.prmtop -t water_solvated-qmmmnvt.netcdf -a 0 -i 1 -s 4``
+``autosolvate clustergen -f water_solvated.prmtop -t water_solvated-qmmmnvt.netcdf -a 0 -i 10 -s 6``
 
 As Autosolvate is running, you will notice this line now includes the list of the 10 frames that the clusters will be extracted from::
 
   extracting from frames: [0, 10, 20, 30, 40, 50, 60, 70, 80, 90]
 
+.. note::
+
+  If you were not able to run the clustergen command above, you can download one of the microsolvated clusters with 6 Å solvent shell size here:
+
+  :download:`water_solvated-cutoutn-0.xyz <_data/water_solvated-cutoutn-0.xyz>`
+
+
 .. warning::
 
    The naming of the microsolvated clusters is based on the name of the .prmtop file, not the trajectory file, so the names will not change between runs. This means that if you run the clustergen command twice, *the new coordinates will overwrite the old ones* (if the frame number is the same). Therefore, if you want to extract clusters from multiple MD steps (like QMMM heat and QMMM NVT), you need to either move or rename the files before you run the command again.