Merge branch 'main' of github.com:alec-glisman/Polyelectrolyte-Surfac…

…e-Adsorption
alec-glisman · Mar 13, 2024 · 91d7c7d · 91d7c7d
2 parents 06afc74 + 1f5d5ce
commit 91d7c7d
Show file tree

Hide file tree

Showing 9 changed files with 334 additions and 33 deletions.
diff --git a/.vscode/settings.json b/.vscode/settings.json
@@ -75,6 +75,7 @@
         "Mdrun",
         "METAD",
         "metadynamics",
+        "mindepth",
         "mpic",
         "mpicc",
         "mpiexec",

diff --git a/scripts/method/equilibration.sh b/scripts/method/equilibration.sh
@@ -95,12 +95,6 @@ else
             sed -i 's/^rcoulomb.*/rcoulomb = 0.7/g' "${sim_name}.mdp" || exit 1
             sed -i 's/^rvdw.*/rvdw = 0.7/g' "${sim_name}.mdp" || exit 1
         fi
-        # add vacuum parameters to mdp file
-        if [[ "${VACUUM_HEIGHT}" -gt 0 ]]; then
-            sed -i 's/^ewald-geometry .*/ewald-geometry            = 3dc/g' "${sim_name}.mdp" || exit 1
-            sed -i 's/^pbc .*/pbc                       = xy/g' "${sim_name}.mdp" || exit 1
-            sed -i 's/^nwall .*/nwall                     = 2/g' "${sim_name}.mdp" || exit 1
-        fi
 
         # make tpr file for NVT equilibration
         "${GMX_BIN}" -nocopyright grompp \
@@ -193,12 +187,6 @@ else
             sed -i 's/^rcoulomb.*/rcoulomb = 0.7/g' "${sim_name}.mdp" || exit 1
             sed -i 's/^rvdw.*/rvdw = 0.7/g' "${sim_name}.mdp" || exit 1
         fi
-        # add vacuum parameters to mdp file
-        if [[ "${VACUUM_HEIGHT}" -gt 0 ]]; then
-            sed -i 's/^ewald-geometry .*/ewald-geometry            = 3dc/g' "${sim_name}.mdp" || exit 1
-            sed -i 's/^pbc .*/pbc                       = xy/g' "${sim_name}.mdp" || exit 1
-            sed -i 's/^nwall .*/nwall                     = 2/g' "${sim_name}.mdp" || exit 1
-        fi
 
         # make tpr file
         "${GMX_BIN}" -nocopyright grompp \
@@ -291,6 +279,24 @@ for i in "${!previous_sim_gro_box_dimensions_array[@]}"; do
     echo "DEBUG: Percent change in dimension ${i}: ${percent_change}"'%'
 done
 
+# ##############################################################################
+# Add Vacuum Layer #############################################################
+# ##############################################################################
+if [[ "${VACUUM_HEIGHT}" -eq 0 ]]; then
+    echo "INFO: Skipping vacuum layer"
+else
+    echo "INFO: Adding vacuum layer"
+    {
+        # get box height from last line of gro file
+        box_dim="$(tail -n 1 "${sim_name}.gro")"
+        z_height="$(echo "${box_dim}" | awk '{print $3}')"
+        z_vacuum_height="$(bc <<<"scale=5; ${z_height} + ${VACUUM_HEIGHT}")"
+
+        # replace z-dimension of box with vacuum layer by string replacement of 3rd column in last line of gro file
+        sed -i "s/${z_height}/${z_vacuum_height}/g" "${sim_name}.gro"
+    } >>"${log_file}" 2>&1
+fi
+
 # #######################################################################################
 # Production equilibration ##############################################################
 # #######################################################################################

diff --git a/scripts/method/initialization.sh b/scripts/method/initialization.sh
@@ -106,14 +106,6 @@ echo "INFO: Copying input files to working directory"
         sed -i 's/^rcoulomb.*/rcoulomb = 0.7/g' "mdin.mdp"
         sed -i 's/^rvdw.*/rvdw = 0.7/g' "mdin.mdp"
     fi
-
-    # add vacuum parameters to mdp file
-    if [[ "${VACUUM_HEIGHT}" -gt 0 ]]; then
-        sed -i 's/^ewald-geometry .*/ewald-geometry            = 3dc/g' "mdin.mdp"
-        sed -i 's/^pbc .*/pbc                       = xy/g' "mdin.mdp"
-        sed -i 's/^nwall .*/nwall                     = 2/g' "mdin.mdp"
-    fi
-
 } >>"${log_file}" 2>&1
 
 # ##############################################################################
@@ -382,17 +374,6 @@ System
 EOF
 } >>"${log_file}" 2>&1
 
-# ##############################################################################
-# Add Vacuum Layer #############################################################
-# ##############################################################################
-echo "INFO: Adding vacuum layer"
-{
-    # increase z-dimension of box with vacuum layer by string replacement of 3rd column in last line of gro file
-    z_box_height="$(bc <<<"scale=5; ${BOX_HEIGHT} * 1.00000")"
-    z_box_vacuum_height="$(bc <<<"scale=5; ${z_box_height} + ${VACUUM_HEIGHT}")"
-    sed -i "s/${z_box_height}/${z_box_vacuum_height}/g" "${sim_name}.gro"
-} >>"${log_file}" 2>&1
-
 # ##############################################################################
 # Make Index File ##############################################################
 # ##############################################################################

diff --git a/submission/1eqbm_7.sh b/submission/1eqbm_7.sh
@@ -0,0 +1,64 @@
+#!/usr/bin/env bash
+# Created by Alec Glisman (GitHub: @alec-glisman) on January 25th, 2023
+
+#SBATCH --time=2-00:00:00
+
+# Slurm: Node configuration
+#SBATCH --partition=all --qos=dow --account=dow
+#SBATCH --nodes=1 --ntasks-per-node=16 --mem=4G
+#SBATCH --gres=gpu:1 --gpu-bind=closest
+
+# Slurm: Runtime I/O
+#SBATCH --mail-user=slurm.notifications@gmail.com --mail-type=BEGIN,END,FAIL
+#SBATCH --output=logs/jobid_%j-node_%N-%x.log --error=logs/jobid_%j-node_%N-%x.log
+
+# REVIEW: find all input files with command:
+# find . -maxdepth 1 -type f -print0 | xargs -0 -I {} echo "'{}'" | sed "s|'./|'|g" | sort
+
+# built-in shell options
+set -o errexit # exit when a command fails
+set -o nounset # exit when script tries to use undeclared variables
+
+# simulation path variables
+proj_base_dir="$(pwd)/.."
+scripts_dir="${proj_base_dir}/scripts"
+params_dir="${proj_base_dir}/submission/input/7-vacuum-systems"
+mapfile -t input_globals < <(find "${params_dir}" -mindepth 1 -maxdepth 1 -type f -print0 | xargs -0 -I {} basename {} | sort)
+
+echo "INFO: Found ${#input_globals[@]} input global directories"
+for idx in "${!input_globals[@]}"; do
+    echo "DEBUG: ${idx}: ${input_globals[idx]}"
+done
+
+# argument parsing
+# 1st argument: input global indices to run
+idxs=()
+if [ $# -eq 0 ]; then
+    mapfile -t idxs < <(seq 0 $((${#input_globals[@]} - 1)))
+else
+    idxs=("$@")
+fi
+
+# check if input indices are valid
+for idx in "${idxs[@]}"; do
+    if [ "${idx}" -lt 0 ] || [ "${idx}" -ge ${#input_globals[@]} ]; then
+        echo "ERROR: invalid input index: ${idx}"
+        exit 1
+    fi
+done
+
+# filter input globals
+input_globals=("${input_globals[@]:${idxs[0]}:${#idxs[@]}}")
+echo "DEBUG: Running input global directories:" "${input_globals[@]}"
+
+# start script
+date_time=$(date +"%Y-%m-%d %T")
+echo "START: ${date_time}"
+
+parallel --link --keep-order --ungroup --halt-on-error '2' --jobs '1' \
+    "${scripts_dir}/run.sh" "${params_dir}/{1}" --equilibrate \
+    ::: "${input_globals[@]}"
+
+# end script
+date_time=$(date +"%Y-%m-%d %T")
+echo "END: ${date_time}"
diff --git a/submission/input/7-vacuum-systems/1.0_104calcite_9nm_crystal_10nm_height.sh b/submission/input/7-vacuum-systems/1.0_104calcite_9nm_crystal_10nm_height.sh
@@ -0,0 +1,83 @@
+#!/usr/bin/env bash
+# -*- coding: utf-8 -*-
+#
+# Author     : Alec Glisman (GitHub: @alec-glisman)
+# Date       : 2023-08-30
+# Description: Script to set global variables and preferences for the simulation.
+# Notes      : Script should only be called from the main run.sh script.
+
+# Hardware ####################################################################
+
+export CPU_THREADS='-1' # number of CPU threads to use (-1 = all available)
+export PIN_OFFSET='-1'  # offset for CPU thread pinning (-1 = no offset)
+export GPU_IDS='-1'     # GPU device(s) to use (0 = first GPU, 01 = first two GPUs)
+
+# System components ###########################################################
+
+# tag for system
+export TAG_JOBID="7.1.0-00idx" # tag to append to system name
+
+# statistical mechanics
+export PRODUCTION_ENSEMBLE='NVT' # {NVT, NPT}
+export TEMPERATURE_K='300'       # temperature in Kelvin
+export PRESSURE_BAR='1'          # pressure in bar
+
+# integration
+export INTEGRATION_NS='500' # [ns] final simulation time for production run {100, 500}
+
+# polyelectrolyte chemistry
+export MONOMER="Acr" # Dominant monomer: {Acr, Acn, Asp, Glu, Ace, Alc}
+export BLOCK=""      # Block copolymer: {iiia, iiiiiiia, iiic, iiiiiic} for {Acr, Acn}
+
+# calcium carbonate crystal surface
+export CRYSTAL="calcite" # {calcite, aragonite, vaterite}
+export SURFACE="104"     # Miller index of crystal surface {104, 001}
+
+# system size
+export SURFACE_SIZE='9'     # size of crystal surface in nm {3, 5, 8, 9, 10, 11, 12, 13}
+export BOX_HEIGHT='10'      # height of simulation box in nm
+export VACUUM_HEIGHT='30'   # height of vacuum layer in nm
+export PDB_BULK_ZMIN='4.62' # z-coordinate of bottom of bulk part of crystal in nm in PDB file
+export PDB_BULK_ZMAX='5.09' # z-coordinate of top of bulk part of crystal in nm in PDB file
+
+# number of each component
+export N_MONOMER='0'   # number of monomers in chain {1, 2, 5, 8, 16, 32}
+export N_CHAIN='0'     # number of chains
+export N_CARBONATE='0' # number of aqueous carbonate ions
+export N_SODIUM='0'    # number of aqueous sodium ions
+export N_CALCIUM='0'   # number of aqueous calcium ions
+export N_CHLORINE='0'  # number of aqueous chlorine ions
+
+# Enhanced sampling ###########################################################
+
+# hamiltonian replica exchange
+export HREMD_N_REPLICA='32' # number of replicas in HREMD simulations
+export HREMD_N_STEPS='1000' # number of steps between replica exchange attempts
+
+# OneOPES replica exchange
+export ONEOPES_N_REPLICA='8'        # number of replicas in OneOPES simulations
+export ONEOPES_N_STEPS='1000'       # number of steps between replica exchange attempts
+export ONEOPES_LARGE_BARRIER='30'   # [kJ/mol] large barrier height for OneOPES replica exchange
+export ONEOPES_SMALL_BARRIER='5'    # [kJ/mol] small barrier height for OneOPES replica exchange
+export ONEOPES_REPLICA_2_TEMP='300' # [K] max OPES MultiTherm temperature of replica 2
+export ONEOPES_REPLICA_3_TEMP='310' # [K] max OPES MultiTherm temperature of replica 3
+export ONEOPES_REPLICA_4_TEMP='325' # [K] max OPES MultiTherm temperature of replica 4
+export ONEOPES_REPLICA_5_TEMP='340' # [K] max OPES MultiTherm temperature of replica 5
+export ONEOPES_REPLICA_6_TEMP='355' # [K] max OPES MultiTherm temperature of replica 6
+export ONEOPES_REPLICA_7_TEMP='370' # [K] max OPES MultiTherm temperature of replica 7
+
+# well-tempered metadynamics
+export METAD_BIASFACTOR='8'       # bias factor for C.V. effective temperature
+export METAD_HEIGHT='1.0'         # [kJ/mol] initial height of Gaussians (kT = 2.48 kJ/mol at 298 K)
+export METAD_SIGMA='0.025'        # width of Gaussians, set to 0.33–0.5 of estimated fluctuation
+export METAD_GRID_SPACING='0.005' # width of bins in the meta-dynamics grid
+export METAD_GRID_MIN='0'         # minimum grid point for Gaussian deposition
+export METAD_GRID_MAX='10'        # maximum grid point for Gaussian deposition
+export METAD_PACE='500'           # [steps] between deposition of Gaussians
+
+# harmonic restraints
+export PE_WALL_MIN='0.3'      # z-coordinate of lower wall in nm
+export PE_WALL_MAX='4.0'      # z-coordinate of upper wall in nm
+export PE_WALL_MAX_EQBM='2.0' # z-coordinate of upper wall in nm during equilibration
+export ATOM_REFERENCE='7310'  # atom number of reference atom for harmonic restraints (1 = first atom)
+export ATOM_OFFSET='-0.305'   # z-coordinate offset of reference atom from crystal surface in nm
diff --git a/submission/input/7-vacuum-systems/2.0_64CaCl2_104calcite_9nm_crystal_10nm_height.sh b/submission/input/7-vacuum-systems/2.0_64CaCl2_104calcite_9nm_crystal_10nm_height.sh
@@ -0,0 +1,83 @@
+#!/usr/bin/env bash
+# -*- coding: utf-8 -*-
+#
+# Author     : Alec Glisman (GitHub: @alec-glisman)
+# Date       : 2023-08-30
+# Description: Script to set global variables and preferences for the simulation.
+# Notes      : Script should only be called from the main run.sh script.
+
+# Hardware ####################################################################
+
+export CPU_THREADS='-1' # number of CPU threads to use (-1 = all available)
+export PIN_OFFSET='-1'  # offset for CPU thread pinning (-1 = no offset)
+export GPU_IDS='-1'     # GPU device(s) to use (0 = first GPU, 01 = first two GPUs)
+
+# System components ###########################################################
+
+# tag for system
+export TAG_JOBID="7.2.0-01idx" # tag to append to system name
+
+# statistical mechanics
+export PRODUCTION_ENSEMBLE='NVT' # {NVT, NPT}
+export TEMPERATURE_K='300'       # temperature in Kelvin
+export PRESSURE_BAR='1'          # pressure in bar
+
+# integration
+export INTEGRATION_NS='500' # [ns] final simulation time for production run {100, 500}
+
+# polyelectrolyte chemistry
+export MONOMER="Acr" # Dominant monomer: {Acr, Acn, Asp, Glu, Ace, Alc}
+export BLOCK=""      # Block copolymer: {iiia, iiiiiiia, iiic, iiiiiic} for {Acr, Acn}
+
+# calcium carbonate crystal surface
+export CRYSTAL="calcite" # {calcite, aragonite, vaterite}
+export SURFACE="104"     # Miller index of crystal surface {104, 001}
+
+# system size
+export SURFACE_SIZE='9'     # size of crystal surface in nm {3, 5, 8, 9, 10, 11, 12, 13}
+export BOX_HEIGHT='10'      # height of simulation box in nm
+export VACUUM_HEIGHT='30'   # height of vacuum layer in nm
+export PDB_BULK_ZMIN='4.62' # z-coordinate of bottom of bulk part of crystal in nm in PDB file
+export PDB_BULK_ZMAX='5.09' # z-coordinate of top of bulk part of crystal in nm in PDB file
+
+# number of each component
+export N_MONOMER='0'    # number of monomers in chain {1, 2, 5, 8, 16, 32}
+export N_CHAIN='0'      # number of chains
+export N_CARBONATE='0'  # number of aqueous carbonate ions
+export N_SODIUM='0'     # number of aqueous sodium ions
+export N_CALCIUM='64'   # number of aqueous calcium ions
+export N_CHLORINE='128' # number of aqueous chlorine ions
+
+# Enhanced sampling ###########################################################
+
+# hamiltonian replica exchange
+export HREMD_N_REPLICA='32' # number of replicas in HREMD simulations
+export HREMD_N_STEPS='1000' # number of steps between replica exchange attempts
+
+# OneOPES replica exchange
+export ONEOPES_N_REPLICA='8'        # number of replicas in OneOPES simulations
+export ONEOPES_N_STEPS='1000'       # number of steps between replica exchange attempts
+export ONEOPES_LARGE_BARRIER='30'   # [kJ/mol] large barrier height for OneOPES replica exchange
+export ONEOPES_SMALL_BARRIER='5'    # [kJ/mol] small barrier height for OneOPES replica exchange
+export ONEOPES_REPLICA_2_TEMP='300' # [K] max OPES MultiTherm temperature of replica 2
+export ONEOPES_REPLICA_3_TEMP='310' # [K] max OPES MultiTherm temperature of replica 3
+export ONEOPES_REPLICA_4_TEMP='325' # [K] max OPES MultiTherm temperature of replica 4
+export ONEOPES_REPLICA_5_TEMP='340' # [K] max OPES MultiTherm temperature of replica 5
+export ONEOPES_REPLICA_6_TEMP='355' # [K] max OPES MultiTherm temperature of replica 6
+export ONEOPES_REPLICA_7_TEMP='370' # [K] max OPES MultiTherm temperature of replica 7
+
+# well-tempered metadynamics
+export METAD_BIASFACTOR='8'       # bias factor for C.V. effective temperature
+export METAD_HEIGHT='1.0'         # [kJ/mol] initial height of Gaussians (kT = 2.48 kJ/mol at 298 K)
+export METAD_SIGMA='0.025'        # width of Gaussians, set to 0.33–0.5 of estimated fluctuation
+export METAD_GRID_SPACING='0.005' # width of bins in the meta-dynamics grid
+export METAD_GRID_MIN='0'         # minimum grid point for Gaussian deposition
+export METAD_GRID_MAX='10'        # maximum grid point for Gaussian deposition
+export METAD_PACE='500'           # [steps] between deposition of Gaussians
+
+# harmonic restraints
+export PE_WALL_MIN='0.0'      # z-coordinate of lower wall in nm
+export PE_WALL_MAX='20.0'     # z-coordinate of upper wall in nm
+export PE_WALL_MAX_EQBM='2.0' # z-coordinate of upper wall in nm during equilibration
+export ATOM_REFERENCE='7310'  # atom number of reference atom for harmonic restraints (1 = first atom)
+export ATOM_OFFSET='-0.305'   # z-coordinate offset of reference atom from crystal surface in nm