diff --git a/Dockerfile.gromacs_amber b/Dockerfile.gromacs_amber new file mode 100644 index 0000000..67c5abd --- /dev/null +++ b/Dockerfile.gromacs_amber @@ -0,0 +1,192 @@ +# syntax=docker/dockerfile:1 +# NOTE: Building this image require's docker version >= 23.0. +# +# For reference: +# - https://docs.docker.com/build/dockerfile/frontend/#stable-channel +ARG TAG_VERSION="12.4.1" +FROM nvidia/cuda:${TAG_VERSION}-cudnn-devel-ubuntu22.04 +ARG HTTP_PROXY +ARG HTTPS_PROXY +ENV http_proxy=${HTTP_PROXY} +ENV https_proxy=${HTTPS_PROXY} +ARG DEBIAN_FRONTEND="noninteractive" +ENV DEBIAN_FRONTEND=${DEBIAN_FRONTEND} +ARG ROOT_PASSWD="root" +ENV ROOT_PASSWD=${ROOT_PASSWD} +ENV SSH_PORT=2222 +WORKDIR /root +SHELL ["/bin/bash", "-c"] + +# base tools +RUN < ~/.ssh/config +cp /etc/ssh/sshd_config /etc/ssh/sshd_config.bak +sed -i 's/#PermitRootLogin prohibit-password/PermitRootLogin yes/' /etc/ssh/sshd_config +sed -i 's/#PasswordAuthentication yes/PasswordAuthentication yes/' /etc/ssh/sshd_config +sed -i 's/#PubkeyAuthentication yes/PubkeyAuthentication yes/' /etc/ssh/sshd_config +sed -i 's/^\(\s*\)GSSAPIAuthentication yes/\1GSSAPIAuthentication no/' /etc/ssh/ssh_config +sed -i "s/^#Port 22/Port ${SSH_PORT}/" /etc/ssh/sshd_config +sudo sed -i "s/# Port 22/Port ${SSH_PORT}/" /etc/ssh/ssh_config +ssh-keygen -t rsa -b 4096 -f /root/.ssh/id_rsa -N "" <<< y +cat ~/.ssh/id_rsa.pub >> ~/.ssh/auth +cat /root/.ssh/id_rsa.pub >> /root/.ssh/authorized_keys +cat /root/.ssh/id_rsa.pub >> /root/.ssh/authorized_keys2 +chmod 600 /root/.ssh/authorized_keys +chmod 600 /root/.ssh/authorized_keys2 +mkdir /var/run/sshd +echo "root:${ROOT_PASSWD}" | chpasswd +mkdir -p ~/.pip +# install pixi +curl -fsSL https://pixi.sh/install.sh | bash +EOT + +ARG FFTW_VERSION="3.3.10" +ENV FFTW_VERSION=${FFTW_VERSION} +ENV PATH=/usr/local/fftw:$PATH +# 安装fftw +RUN < ./test_mpi_cuda.cu +#include +#include +#include + +__global__ void hello_cuda() { + printf("Hello from CUDA kernel! Thread id: %d\n", threadIdx.x); +} + +int main(int argc, char **argv) { + MPI_Init(&argc, &argv); + + int rank; + MPI_Comm_rank(MPI_COMM_WORLD, &rank); + + printf("Hello from MPI process %d!\n", rank); + + // Launch CUDA kernel + hello_cuda<<<1, 10>>>(); + cudaDeviceSynchronize(); // Wait for the CUDA kernel to finish + + MPI_Finalize(); + return 0; +} +EOF +nvcc -o test_mpi_cuda test_mpi_cuda.cu -I${CUDA_HOME}/include -I${MPI_HOME}/include -L${MPI_HOME}/lib -lcudart -lmpi +# mpirun --allow-run-as-root -np 2 ./test_mpi_cuda +EOT + +# 安装plumed +ARG PLUMED_VERSION="2.9.1" +ENV PLUMED_VERSION=${PLUMED_VERSION} +ENV LD_LIBRARY_PATH=/usr/local/plumed/lib:$LD_LIBRARY_PATH +ENV PATH=/usr/local/plumed:/usr/local/plumed/bin:$PATH +RUN <> /root/.bashrc +EOT + +COPY file/Amber24.tar.bz2 file/AmberTools24.tar.bz2 /root + +ENV DOWNLOAD_MINICONDA="False" +# install ambertools +RUN < ions.mdp +; ions.mdp - used as input into grompp to generate ions.tpr +; Parameters describing what to do, when to stop and what to save +integrator = steep ; Algorithm （steep = steepest descent minimization） +emtol = 1000.0 ; Stop minimization when the maximum force < 1000.0 kJ/mol/nm +emstep = 0.01 ; Minimization step size +nsteps = 50000 ; Maximum number of （minimization） steps to perform +; Parameters describing how to find the neighbors of each atom and how to calculate the interactions +nstlist = 1 ; Frequency to update the neighbor list and long range forces +cutoff-scheme = Verlet ; Buffered neighbor searching +ns_type = grid ; Method to determine neighbor list （simple, grid） +coulombtype = cutoff ; Treatment of long range electrostatic interactions +rcoulomb = 1.0 ; Short-range electrostatic cut-off +rvdw = 1.0 ; Short-range Van der Waals cut-off +pbc = xyz ; Periodic Boundary Conditions in all 3 dimensions +EOF +gmx_mpi grompp -f ions.mdp -c $NAME-solv.gro -p topol.top -o ions.tpr -maxwarn 1 +echo SOL | gmx_mpi genion -s ions.tpr -o $NAME-solv-ions.gro -p topol.top -pname NA -nname CL -conc 0.125 -neutral +# energy minimization of the structure in solvate # ! minim.mdp add by manual +# --- minim.mdp file content --- # +cat << EOF > minim.mdp +; minim.mdp - used as input into grompp to generate em.tpr +; Parameters describing what to do, when to stop and what to save +integrator = steep ; Algorithm （steep = steepest descent minimization） +emtol = 1000.0 ; Stop minimization when the maximum force < 1000.0 kJ/mol/nm +emstep = 0.01 ; Minimization step size +nsteps = 50000 ; Maximum number of （minimization） steps to perform +; Parameters describing how to find the neighbors of each atom and how to calculate the interactions +nstlist = 1 ; Frequency to update the neighbor list and long range forces +cutoff-scheme = Verlet ; Buffered neighbor searching +ns_type = grid ; Method to determine neighbor list （simple, grid） +coulombtype = PME ; Treatment of long range electrostatic interactions +rcoulomb = 1.0 ; Short-range electrostatic cut-off +rvdw = 1.0 ; Short-range Van der Waals cut-off +pbc = xyz ; Periodic Boundary Conditions in all 3 dimensions +EOF +gmx_mpi grompp -f minim.mdp -c $NAME-solv-ions.gro -p topol.top -o em.tpr +gmx_mpi mdrun -v -deffnm em +# optional em, you will need the Xmgrace plotting too +#gmx_mpi energy -f em.edr -o potential.xvg +#position restrain +# gmx_mpi grompp -f posre.mdp -c em.gro -p topol.top -o posre.tpr -r em.gro +# gmx_mpi mdrun -v -deffnm posre +# nvt +# gmx_mpi grompp -f nvt.mdp -c posre.gro -p topol.top -o nvt.tpr +# --- nvt.mdp file content --- # +cat << EOF > nvt.mdp +title = OPLS Lysozyme NVT equilibration +define = -DPOSRES ; position restrain the protein +; Run parameters +integrator = md ; leap-frog integrator +nsteps = 50000 ; 2 * 50000 = 100 ps +dt = 0.002 ; 2 fs +; Output control +nstxout = 500 ; save coordinates every 1.0 ps +nstvout = 500 ; save velocities every 1.0 ps +nstenergy = 500 ; save energies every 1.0 ps +nstlog = 500 ; update log file every 1.0 ps +; Bond parameters +continuation = no ; first dynamics run +constraint_algorithm = lincs ; holonomic constraints +constraints = h-bonds ; bonds involving H are constrained +lincs_iter = 1 ; accuracy of LINCS +lincs_order = 4 ; also related to accuracy +; Nonbonded settings +cutoff-scheme = Verlet ; Buffered neighbor searching +ns_type = grid ; search neighboring grid cells +nstlist = 10 ; 20 fs, largely irrelevant with Verlet +rcoulomb = 1.0 ; short-range electrostatic cutoff （in nm） +rvdw = 1.0 ; short-range van der Waals cutoff （in nm） +DispCorr = EnerPres ; account for cut-off vdW scheme +; Electrostatics +coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics +pme_order = 4 ; cubic interpolation +fourierspacing = 0.16 ; grid spacing for FFT +; Temperature coupling is on +tcoupl = V-rescale ; modified Berendsen thermostat +tc-grps = Protein Non-Protein ; two coupling groups - more accurate +tau_t = 0.1 0.1 ; time constant, in ps +ref_t = 300 300 ; reference temperature, one for each group, in K +; Pressure coupling is off +pcoupl = no ; no pressure coupling in NVT +; Periodic boundary conditions +pbc = xyz ; 3-D PBC +; Velocity generation +gen_vel = yes ; assign velocities from Maxwell distribution +gen_temp = 300 ; temperature for Maxwell distribution +gen_seed = -1 ; generate a random seed +EOF +gmx_mpi grompp -f nvt.mdp -c em.gro -r em.gro -p topol.top -o nvt.tpr +gmx_mpi mdrun -v -deffnm nvt +# optional : Let's analyze the temperature progression, again using energy: +# gmx_mpi energy -f nvt.edr -o temperature.xvg +# npt +# gmx_mpi grompp -f npt.mdp -c nvt.gro -t nvt.cpt -p topol.top -o npt.tpr +# --- npt.mdp file content --- # +cat << EOF > npt.mdp +title = OPLS Lysozyme NPT equilibration +define = -DPOSRES ; position restrain the protein +; Run parameters +integrator = md ; leap-frog integrator +nsteps = 50000 ; 2 * 50000 = 100 ps +dt = 0.002 ; 2 fs +; Output control +nstxout = 500 ; save coordinates every 1.0 ps +nstvout = 500 ; save velocities every 1.0 ps +nstenergy = 500 ; save energies every 1.0 ps +nstlog = 500 ; update log file every 1.0 ps +; Bond parameters +continuation = yes ; Restarting after NVT +constraint_algorithm = lincs ; holonomic constraints +constraints = h-bonds ; bonds involving H are constrained +lincs_iter = 1 ; accuracy of LINCS +lincs_order = 4 ; also related to accuracy +; Nonbonded settings +cutoff-scheme = Verlet ; Buffered neighbor searching +ns_type = grid ; search neighboring grid cells +nstlist = 10 ; 20 fs, largely irrelevant with Verlet scheme +rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm) +rvdw = 1.0 ; short-range van der Waals cutoff (in nm) +DispCorr = EnerPres ; account for cut-off vdW scheme +; Electrostatics +coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics +pme_order = 4 ; cubic interpolation +fourierspacing = 0.16 ; grid spacing for FFT +; Temperature coupling is on +tcoupl = V-rescale ; modified Berendsen thermostat +tc-grps = Protein Non-Protein ; two coupling groups - more accurate +tau_t = 0.1 0.1 ; time constant, in ps +ref_t = 300 300 ; reference temperature, one for each group, in K +; Pressure coupling is on +pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT +pcoupltype = isotropic ; uniform scaling of box vectors +tau_p = 2.0 ; time constant, in ps +ref_p = 1.0 ; reference pressure, in bar +compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 +refcoord_scaling = com +; Periodic boundary conditions +pbc = xyz ; 3-D PBC +; Velocity generation +gen_vel = no ; Velocity generation is off +EOF +gmx_mpi grompp -f npt.mdp -c nvt.gro -r nvt.gro -t nvt.cpt -p topol.top -o npt.tpr +gmx_mpi mdrun -v -deffnm npt +# Optional: Let's analyze the pressure progression, again using energy: type 18 0 +# gmx energy -f npt.edr -o pressure.xvg +# Optional: Let's take a look at density as well, this time using energy and entering "24 0" at the prompt. +# gmx energy -f npt.edr -o density.xvg +# md +# --- md.mdp file content --- # +cat << EOF > ${MDRUN_NAME}.mdp +title = OPLS Lysozyme NPT equilibration +; Run parameters +integrator = md ; leap-frog integrator +nsteps = ${NSTEPS} ; steps for simulation +dt = ${DT} ; time step in fs +; Output control +nstxout = 0 ; suppress bulky .trr file by specifying +nstvout = 0 ; 0 for output frequency of nstxout, +nstfout = 0 ; nstvout, and nstfout +nstenergy = 5000 ; save energies every 10.0 ps +nstlog = 5000 ; update log file every 10.0 ps +nstxout-compressed = 5000 ; save compressed coordinates every 10.0 ps +compressed-x-grps = System ; save the whole system +; Bond parameters +continuation = yes ; Restarting after NPT +constraint_algorithm = lincs ; holonomic constraints +constraints = h-bonds ; bonds involving H are constrained +lincs_iter = 1 ; accuracy of LINCS +lincs_order = 4 ; also related to accuracy +; Neighborsearching +cutoff-scheme = Verlet ; Buffered neighbor searching +ns_type = grid ; search neighboring grid cells +nstlist = 10 ; 20 fs, largely irrelevant with Verlet scheme +rcoulomb = 1.0 ; short-range electrostatic cutoff （in nm） +rvdw = 1.0 ; short-range van der Waals cutoff （in nm） +; Electrostatics +coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics +pme_order = 4 ; cubic interpolation +fourierspacing = 0.16 ; grid spacing for FFT +; Temperature coupling is on +tcoupl = V-rescale ; modified Berendsen thermostat +tc-grps = Protein Non-Protein ; two coupling groups - more accurate +tau_t = 0.1 0.1 ; time constant, in ps +ref_t = 300 300 ; reference temperature, one for each group, in K +; Pressure coupling is on +pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT +pcoupltype = isotropic ; uniform scaling of box vectors +tau_p = 2.0 ; time constant, in ps +ref_p = 1.0 ; reference pressure, in bar +compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 +; Periodic boundary conditions +pbc = xyz ; 3-D PBC +; Dispersion correction +DispCorr = EnerPres ; account for cut-off vdW scheme +; Velocity generation +gen_vel = no ; Velocity generation is off +EOF + +# Generate GROMACS .tpr file for the simulation +gmx_mpi grompp -f ${MDRUN_NAME}.mdp -c npt.gro -t npt.cpt -p topol.top -o ${TPR_FILE} + +# Run the simulation +gmx_mpi mdrun -deffnm ${MDRUN_NAME} + +# extra ndx file , select protein +echo -e "1\nq" | gmx_mpi make_ndx -f ${MDRUN_NAME}.gro -o ${NDX_FILE} +# echo -e "1\nq" | gmx_mpi make_ndx -f md.gro -o index.ndx + +# Create extraction output directory +mkdir -p ${OUTPUT_FOLDER} + +# Create temp output directory +mkdir -p ${TEMP_FOLDER} + +echo -e "1\nq" | gmx_mpi trjconv -dt ${EXTRACT_EVERY_PS} -s ${TPR_FILE} -f ${XTC_FILE} -n ${NDX_FILE} -pbc mol -o ${TEMP_FOLDER}/temp.xtc + +# echo -e "1\nq" | gmx_mpi trjconv -dt 100 -s md.tpr -f md.xtc -n index.ndx -pbc mol -o temp/temp.xtc + +echo -e "1\n1\n1" | gmx_mpi trjconv -s ${TPR_FILE} -f ${TEMP_FOLDER}/temp.xtc -n ${NDX_FILE} -center -fit rot+trans -o ${TEMP_FOLDER}/traj_show.xtc + +# echo -e "1\n1\n1" | gmx_mpi trjconv -s md.tpr -f temp/temp.xtc -n index.ndx -center -fit rot+trans -o temp/traj_show.xtc + +echo -e "1\n1\n1" | gmx_mpi trjconv -s ${TPR_FILE} -f ${TEMP_FOLDER}/temp.xtc -n ${NDX_FILE} -center -fit rot+trans -b 0 -e 0 -o ${TEMP_FOLDER}/tarj_show.pdb + +# echo -e "1\n1\n1" | gmx_mpi trjconv -s md.tpr -f temp/temp.xtc -n index.ndx -center -fit rot+trans -b 0 -e 0 -o temp/tarj_show.pdb + +# Group 1 ( Protein) + +# --- +# Step 1: Extract frames every 1000 ps +gmx_mpi trjconv -s ${TPR_FILE} -f ${XTC_FILE} -o ${OUTPUT_FOLDER}/${NO_PBC_XTC_FILE} -dt ${EXTRACT_EVERY_PS} -pbc mol <