Simulation and analysis scripts for cubic unit cells workflow. Note: all details of implemented scan_space.py and build_hull.py are not incuded such that only the most important components are shown, but, I am happy to share the full scripts if you reach out to statt@illinois.edu or jamesft2@illinois.edu.

Example workflow:

1. Implement grid or active learning based search in scan_space.py a. Scan_space.py submits jobs to cluster containing thermodynamic_integration.py and all parameters b. Thermodynamic_integration.py runs all requisite thermo integration scripts
2. Run curv.py, update_nrgs_npys.py to compute curvatures, energies, topology measures.
3. Run build_hull.py to build convex hull/phase space. Visualize. Repeat step 1-3.

• Simulations are conducted in HOOMD-blue version 4.8.?
• This relies on our implementation of Martini3 in HOOMD-blue: https://github.com/stattlab/martini3. This package should exist in one's python path to use.

Configured to submit 20 simulations to a cluster for each of the symmetries (calls thermodynamic_integration.py)

Master script which runs all TI scripts sequentially:

• init.py
• equil.py
• add_external_field.py (leg 1)
• add_interactions.py (leg 2)
• remove_external field.py (leg 3)
• stage_4.py (brief stability check) Handles deletion of files in the case of partial running and failure

Initializes unit cell with ghost particles along Im3m Pn3m Ia3d hexagonal or bilayer surfaces... sets aspect ratio of bilayer/hexagonal simulation.

Slowly increase the magnitude of forces from dpd type -> very weak LJ to weak LJ and randomize configurations. All interactions are of the same tiype. SImulate for 360,000 steps.

Over the range of 1,000,000 timesteps, increase the reaction coordinate (alpha) every 10,000 timesteps. Strength of external field (gauss potential) increases linearly with reaction coordinate.Record internal energy of gauss potential.

Increase the strength of the internal interactions over time, increasing the reaction coordinate every 10,000 timesteps. Also decrease the strength of the weak leonard jones potential each timestep. $\epsilon$ = $\epsilon_{LJ}*\eta$. Record internal energy of Martini and non-specific LJ interaction

Render image of the last frame, compute the free energy over all frames, and update an npy which contrains info about the toplogy

Requires:

• Ovito
• scipy
• gsd.hoomd

Assumes

• data is saved to ../data/{hash_value}

Inputs:

• Configured such that when python curv.py is executed, the analysis will be completed for all files within '../data' that end with 'stage_3.gsd'

Outputs:

• Updated .npy files for each sym which contain information about the points sampled {sym}_points.npy, free energies {sym}_values.npy, and final topology {sym}_scores.py

Compute per-lipid area per lipid (local and global), mean curvature, gaussian curvature (local and global) and save these quantities back into a gsd.

Requires:

• Ovito
• networkx
• sklearn
• gsd.hoomd
• freud
• heapq

Assumes:

• Lipid is 6 beads long
• Data is saved to '../data'

Inputs:

• Configured such that when python curv.py is executed, the analysis will be completed for all files within '../data' that end with 'stage_3.gsd'

Outputs:

• .gsd with per-lipid quantities
• text file with formatted average quantities per frame