From 84bd6e245b9b6a93b0b0d597c4a4eab8a237b0ab Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Fri, 26 Aug 2022 10:56:52 +0200 Subject: [PATCH 1/9] including correction files --- transformato/bin/perform_correction.py | 76 ++++++++++++++++++++++ transformato/bin/submit_switching_slurm.sh | 17 +++++ 2 files changed, 93 insertions(+) create mode 100644 transformato/bin/perform_correction.py create mode 100644 transformato/bin/submit_switching_slurm.sh diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py new file mode 100644 index 00000000..64d6671c --- /dev/null +++ b/transformato/bin/perform_correction.py @@ -0,0 +1,76 @@ +# general imports +from endstate_correction.system import create_charmm_system, read_box +from openmm.app import ( + PME, + CharmmParameterSet, + CharmmPsfFile, + PDBFile, +) +from endstate_correction.analysis import plot_endstate_correction_results +import endstate_correction +from endstate_correction.protocoll import perform_endstate_correction, Protocoll +import mdtraj +from openmm import unit + +######################################################## +######################################################## +# ------------ set up the waterbox system -------------- + + +system_name = "ethane" +env = "vacuum" +# define the output directory +output_base = f"." +parameter_base = f"/site/raid3/johannes/free_solv_test/data" +# load the charmm specific files (psf, pdb, rtf, prm and str files) +psf_file = f"../methanol/intst1/lig_in_{env}.psf" +psf = CharmmPsfFile(psf_file) +pdb = PDBFile(f"../methanol/intst1/lig_in_{env}.pdb") +params = CharmmParameterSet( + f"{parameter_base}/{system_name}/waterbox/unk/unk.str", + f"{parameter_base}/{system_name}/waterbox/toppar/top_all36_cgenff.rtf", + f"{parameter_base}/{system_name}/waterbox/toppar/par_all36_cgenff.prm", + f"{parameter_base}/{system_name}/waterbox/toppar/toppar_water_ions.str", +) +# set up the treatment of the system for the specific environment +if env == "waterbox": + psf = read_box(psf, f"{parameter_base}/{system_name}/waterbox/openmm/sysinfo.dat") + +# define region that should be treated with the qml +chains = list(psf.topology.chains()) +ml_atoms = [atom.index for atom in chains[0].atoms()] +# define system +sim = create_charmm_system(psf=psf, parameters=params, env=env, ml_atoms=ml_atoms) + +######################################################## +######################################################## +# ------------------- load samples ---------------------# +n_samples = 500 +n_steps_per_sample = 1250000 +traj_base = f"../methanol/intst1/" +mm_samples = [] +traj = mdtraj.load_dcd( + f"{traj_base}/run_1/lig_in_{env}.dcd", + top=psf_file, +) +if env == "waterbox": + traj.image_molecules() + +mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! + + +#################################################### +# ----------------------- FEP ---------------------- +#################################################### + +fep_protocoll = Protocoll( + method="NEQ", + direction="unidirectional", + sim=sim, + trajectories=[mm_samples], + nr_of_switches=500, # 500 + neq_switching_length=5000, # 5000 +) + +r = perform_endstate_correction(fep_protocoll) +plot_endstate_correction_results(system_name, r, "results_neq_unidirectional.png") diff --git a/transformato/bin/submit_switching_slurm.sh b/transformato/bin/submit_switching_slurm.sh new file mode 100644 index 00000000..53fb4a94 --- /dev/null +++ b/transformato/bin/submit_switching_slurm.sh @@ -0,0 +1,17 @@ +#!/bin/bash +#SBATCH -p lgpu +#SBATCH --gres=gpu +#SBATCH --job-name=switching # Job name +#SBATCH --ntasks=1 # Run on a single CPU +#SBATCH --output=switching_%j.log # Standard output and error log + + +source /home/johannes/miniconda3/etc/profile.d/conda.sh +conda activate endstate + +export OPENMM_PRECISION='mixed' + +pwd; hostname; date + + +python perform_correction.py From 3c29ed10d0d4429ba71d6abc93d141a30d1ee340 Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Mon, 29 Aug 2022 17:38:01 +0200 Subject: [PATCH 2/9] creating option in perfrom_mutation function --- transformato/bin/perform_correction.py | 131 +++++++++++++++---------- transformato/mutate.py | 7 ++ transformato/state.py | 19 ++++ transformato/tests/test_absolute.py | 28 ++++++ 4 files changed, 131 insertions(+), 54 deletions(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index 64d6671c..7538945c 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -1,76 +1,99 @@ # general imports +from os import system from endstate_correction.system import create_charmm_system, read_box from openmm.app import ( PME, CharmmParameterSet, CharmmPsfFile, PDBFile, + CharmmCrdFile, ) from endstate_correction.analysis import plot_endstate_correction_results import endstate_correction from endstate_correction.protocoll import perform_endstate_correction, Protocoll import mdtraj from openmm import unit +import sys ######################################################## ######################################################## # ------------ set up the waterbox system -------------- +system_name = sys.argv[1] -system_name = "ethane" -env = "vacuum" -# define the output directory -output_base = f"." -parameter_base = f"/site/raid3/johannes/free_solv_test/data" -# load the charmm specific files (psf, pdb, rtf, prm and str files) -psf_file = f"../methanol/intst1/lig_in_{env}.psf" -psf = CharmmPsfFile(psf_file) -pdb = PDBFile(f"../methanol/intst1/lig_in_{env}.pdb") -params = CharmmParameterSet( - f"{parameter_base}/{system_name}/waterbox/unk/unk.str", - f"{parameter_base}/{system_name}/waterbox/toppar/top_all36_cgenff.rtf", - f"{parameter_base}/{system_name}/waterbox/toppar/par_all36_cgenff.prm", - f"{parameter_base}/{system_name}/waterbox/toppar/toppar_water_ions.str", -) -# set up the treatment of the system for the specific environment -if env == "waterbox": - psf = read_box(psf, f"{parameter_base}/{system_name}/waterbox/openmm/sysinfo.dat") +def gen_box(psf, crd): + coords = crd.positions -# define region that should be treated with the qml -chains = list(psf.topology.chains()) -ml_atoms = [atom.index for atom in chains[0].atoms()] -# define system -sim = create_charmm_system(psf=psf, parameters=params, env=env, ml_atoms=ml_atoms) + min_crds = [coords[0][0], coords[0][1], coords[0][2]] + max_crds = [coords[0][0], coords[0][1], coords[0][2]] -######################################################## -######################################################## -# ------------------- load samples ---------------------# -n_samples = 500 -n_steps_per_sample = 1250000 -traj_base = f"../methanol/intst1/" -mm_samples = [] -traj = mdtraj.load_dcd( - f"{traj_base}/run_1/lig_in_{env}.dcd", - top=psf_file, -) -if env == "waterbox": - traj.image_molecules() - -mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! - - -#################################################### -# ----------------------- FEP ---------------------- -#################################################### - -fep_protocoll = Protocoll( - method="NEQ", - direction="unidirectional", - sim=sim, - trajectories=[mm_samples], - nr_of_switches=500, # 500 - neq_switching_length=5000, # 5000 -) + for coord in coords: + min_crds[0] = min(min_crds[0], coord[0]) + min_crds[1] = min(min_crds[1], coord[1]) + min_crds[2] = min(min_crds[2], coord[2]) + max_crds[0] = max(max_crds[0], coord[0]) + max_crds[1] = max(max_crds[1], coord[1]) + max_crds[2] = max(max_crds[2], coord[2]) + + boxlx = max_crds[0]-min_crds[0] + boxly = max_crds[1]-min_crds[1] + boxlz = max_crds[2]-min_crds[2] + + psf.setBox(boxlx, boxly, boxlz) + return psf + + +for env in ["waterbox","vacuum"]: + + parameter_base = f"/site/raid3/johannes/free_solv_test/data" + # load the charmm specific files (psf, pdb, rtf, prm and str files) + psf_file = f"../{system_name}/intst1/lig_in_{env}.psf" + psf = CharmmPsfFile(psf_file) + pdb = PDBFile(f"../{system_name}/intst1/lig_in_{env}.pdb") + crd = CharmmCrdFile(f"../{system_name}/intst1/lig_in_{env}.crd") + params = CharmmParameterSet( + f"../{system_name}/intst1/waterbox/unk/unk.str", + f"../toppar/top_all36_cgenff.rtf", + f"../toppar/par_all36_cgenff.prm", + f"../toppar/toppar_water_ions.str", + ) + # set up the treatment of the system for the specific environment + if env == "waterbox": + psf = gen_box(psf, crd) + + # define region that should be treated with the qml + chains = list(psf.topology.chains()) + ml_atoms = [atom.index for atom in chains[0].atoms()] + # define system + sim = create_charmm_system(psf=psf, parameters=params, env=env, ml_atoms=ml_atoms) + + ######################################################## + ######################################################## + # ------------------- load samples ---------------------# + + mm_samples = [] + traj = mdtraj.load_dcd( + f"../{system_name}/intst1/run_1/lig_in_{env}.dcd", + top=psf_file, + ) + if env == "waterbox": + traj.image_molecules() + + mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! + + + #################################################### + # ----------------------- FEP ---------------------- + #################################################### + + fep_protocoll = Protocoll( + method="NEQ", + direction="unidirectional", + sim=sim, + trajectories=[mm_samples], + nr_of_switches=500, # 500 + neq_switching_length=5000, # 5000 + ) -r = perform_endstate_correction(fep_protocoll) -plot_endstate_correction_results(system_name, r, "results_neq_unidirectional.png") + r = perform_endstate_correction(fep_protocoll) + plot_endstate_correction_results(system_name, r, "results_neq_unidirectional.png") diff --git a/transformato/mutate.py b/transformato/mutate.py index 1bf7a57a..14fb4b0b 100644 --- a/transformato/mutate.py +++ b/transformato/mutate.py @@ -18,6 +18,7 @@ from transformato.system import SystemStructure + logger = logging.getLogger(__name__) @@ -70,6 +71,7 @@ def perform_mutations( nr_of_mutation_steps_lj_of_hydrogens: int = 1, nr_of_mutation_steps_lj_of_heavy_atoms: int = 1, nr_of_mutation_steps_cc: int = 5, + endstate_correction: bool = False, ): """Performs the mutations necessary to mutate the physical endstate to the defined common core. @@ -264,6 +266,11 @@ def perform_mutations( mutation=mutation_list["transform"], ) + if endstate_correction: + print("HI sind hier") + i.endstate_correction() + + @dataclass class DummyRegion: diff --git a/transformato/state.py b/transformato/state.py index 1db3735d..aa2cc880 100644 --- a/transformato/state.py +++ b/transformato/state.py @@ -1,6 +1,7 @@ import logging import os import shutil +import glob from io import StringIO from typing import List @@ -43,6 +44,24 @@ def __init__( self.current_step = 1 self.consecutive_runs = consecutive_runs + def endstate_correction(self): + + try: + os.makedirs(f"{self.path}/../endstate_correction") + except: + logger.info(f"Folder for endstate correction exist already") + + # copyt perform_correction.py and the corresponding submit script to the newly created folder + endstate_correction_script_source = f"{self.configuration['bin_dir']}/perform_correction.py" + endstate_correction_script_target = f"{self.path}/../endstate_correction/perform_correction.py" + shutil.copyfile(endstate_correction_script_source, endstate_correction_script_target) + + endstate_correction_script_source = f"{self.configuration['bin_dir']}/submit_switching_slurm.sh" + endstate_correction_script_target = f"{self.path}/../endstate_correction/submit_switching_slurm.sh" + shutil.copyfile(endstate_correction_script_source, endstate_correction_script_target) + + self.configuration["simulation"] + def write_state( self, mutation_conf: List, diff --git a/transformato/tests/test_absolute.py b/transformato/tests/test_absolute.py index 5c0e4095..b9971435 100644 --- a/transformato/tests/test_absolute.py +++ b/transformato/tests/test_absolute.py @@ -112,3 +112,31 @@ def test_compare_mda_and_mdtraj(): mda_results = analyse_asfe_with_module(module="mda") mdtraj_results = analyse_asfe_with_module(module="mdtraj") assert np.isclose(np.average(mda_results), np.average(mdtraj_results)) + + + + +def test_perform_enstate_correction_asfe_system(): + + configuration = load_config_yaml( + config=f"{get_testsystems_dir()}/config/methanol-asfe.yaml", + input_dir=get_testsystems_dir(), + output_dir=get_test_output_dir(), + ) + + s1, mutation_list = create_asfe_system(configuration) + + i = IntermediateStateFactory(system=s1, consecutive_runs=3, configuration=configuration) + + perform_mutations( + configuration=configuration, + nr_of_mutation_steps_charge=2, + i=i, + mutation_list=mutation_list, + endstate_correction=True, + ) + + + assert len(i.output_files) == 7 + assert len((mutation_list)["charge"][0].atoms_to_be_mutated) == 6 + From 8166e320255930db198a052ac34f4167c633a506 Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Tue, 30 Aug 2022 13:25:12 +0200 Subject: [PATCH 3/9] adapting the perform_correction.py script --- transformato/bin/perform_correction.py | 70 +++++++++++++++----------- transformato/mutate.py | 1 - transformato/state.py | 29 +++++++---- 3 files changed, 60 insertions(+), 40 deletions(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index 7538945c..705ba576 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -1,25 +1,23 @@ # general imports -from os import system +import os from endstate_correction.system import create_charmm_system, read_box from openmm.app import ( - PME, CharmmParameterSet, CharmmPsfFile, PDBFile, CharmmCrdFile, ) from endstate_correction.analysis import plot_endstate_correction_results -import endstate_correction -from endstate_correction.protocoll import perform_endstate_correction, Protocoll +from endstate_correction.protocol import perform_endstate_correction, Protocol import mdtraj from openmm import unit -import sys -######################################################## -######################################################## -# ------------ set up the waterbox system -------------- -system_name = sys.argv[1] + +# Variables will be generated by transformato +system_name = NAMEofSYSTEM +tcl = TLC + def gen_box(psf, crd): coords = crd.positions @@ -35,28 +33,38 @@ def gen_box(psf, crd): max_crds[1] = max(max_crds[1], coord[1]) max_crds[2] = max(max_crds[2], coord[2]) - boxlx = max_crds[0]-min_crds[0] - boxly = max_crds[1]-min_crds[1] - boxlz = max_crds[2]-min_crds[2] + boxlx = max_crds[0] - min_crds[0] + boxly = max_crds[1] - min_crds[1] + boxlz = max_crds[2] - min_crds[2] psf.setBox(boxlx, boxly, boxlz) return psf -for env in ["waterbox","vacuum"]: +for env in ["waterbox", "vacuum"]: parameter_base = f"/site/raid3/johannes/free_solv_test/data" - # load the charmm specific files (psf, pdb, rtf, prm and str files) + # load the charmm specific files (psf, pdb, rtf, crd files) psf_file = f"../{system_name}/intst1/lig_in_{env}.psf" psf = CharmmPsfFile(psf_file) pdb = PDBFile(f"../{system_name}/intst1/lig_in_{env}.pdb") crd = CharmmCrdFile(f"../{system_name}/intst1/lig_in_{env}.crd") - params = CharmmParameterSet( - f"../{system_name}/intst1/waterbox/unk/unk.str", - f"../toppar/top_all36_cgenff.rtf", - f"../toppar/par_all36_cgenff.prm", - f"../toppar/toppar_water_ions.str", - ) + + # load forcefiled files (ligand.str and toppar files) + parms = () + file = f"../{system_name}/intst1/{tcl}" + if os.path.isfile(f"{file}.str"): + parms += (f"{file.lower()}.str",) + else: + parms += (f"{file.lower()}_g.rtf",) + parms += (f"{file.lower()}.prm",) + + parms += (f"../toppar/top_all36_cgenff.rtf",) + parms += (f"../toppar/par_all36_cgenff.prm",) + parms += (f"../toppar/toppar_water_ions.str",) + + params = CharmmParameterSet(*parms) + # set up the treatment of the system for the specific environment if env == "waterbox": psf = gen_box(psf, crd) @@ -72,28 +80,30 @@ def gen_box(psf, crd): # ------------------- load samples ---------------------# mm_samples = [] - traj = mdtraj.load_dcd( - f"../{system_name}/intst1/run_1/lig_in_{env}.dcd", - top=psf_file, - ) + if os.path.isfile(f"../{system_name}/intst1/run_1/lig_in_{env}.dcd"): + path_mm = f"../{system_name}/intst1/run_1/lig_in_{env}.dcd" + else: + path_mm = f"../{system_name}/intst1/lig_in_{env}.dcd" + + traj = mdtraj.load_dcd(path_mm, top=psf_file) + if env == "waterbox": traj.image_molecules() - - mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! + mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! #################################################### # ----------------------- FEP ---------------------- #################################################### - fep_protocoll = Protocoll( + fep_protocoll = Protocol( method="NEQ", direction="unidirectional", sim=sim, trajectories=[mm_samples], - nr_of_switches=500, # 500 - neq_switching_length=5000, # 5000 + nr_of_switches=500, # 500 + neq_switching_length=5000, # 5000 ) r = perform_endstate_correction(fep_protocoll) - plot_endstate_correction_results(system_name, r, "results_neq_unidirectional.png") + plot_endstate_correction_results(system_name, r, f"results_neq_unidirectional_{env}.png") diff --git a/transformato/mutate.py b/transformato/mutate.py index 14fb4b0b..23a2e3d4 100644 --- a/transformato/mutate.py +++ b/transformato/mutate.py @@ -267,7 +267,6 @@ def perform_mutations( ) if endstate_correction: - print("HI sind hier") i.endstate_correction() diff --git a/transformato/state.py b/transformato/state.py index aa2cc880..cd52e964 100644 --- a/transformato/state.py +++ b/transformato/state.py @@ -45,22 +45,33 @@ def __init__( self.consecutive_runs = consecutive_runs def endstate_correction(self): - + + logger.info(f"Will create script for endstate correction") try: os.makedirs(f"{self.path}/../endstate_correction") except: logger.info(f"Folder for endstate correction exist already") - # copyt perform_correction.py and the corresponding submit script to the newly created folder + # copy submit script to the newly created folder + submit_switching_script_source = f"{self.configuration['bin_dir']}/submit_switching_slurm.sh" + submit_switchting_script_target = f"{self.path}/../endstate_correction/submit_switching_slurm.sh" + shutil.copyfile(submit_switching_script_source, submit_switchting_script_target) + + # modify the perform_correction file from transformato bin and save it in the endstate_correcition folder endstate_correction_script_source = f"{self.configuration['bin_dir']}/perform_correction.py" endstate_correction_script_target = f"{self.path}/../endstate_correction/perform_correction.py" - shutil.copyfile(endstate_correction_script_source, endstate_correction_script_target) - - endstate_correction_script_source = f"{self.configuration['bin_dir']}/submit_switching_slurm.sh" - endstate_correction_script_target = f"{self.path}/../endstate_correction/submit_switching_slurm.sh" - shutil.copyfile(endstate_correction_script_source, endstate_correction_script_target) - - self.configuration["simulation"] + fin = open(endstate_correction_script_source,"rt") + fout = open(endstate_correction_script_target,"wt") + + for line in fin: + if "NAMEofSYSTEM" in line: + fout.write(line.replace("NAMEofSYSTEM",f'"{self.configuration["system"]["structure1"]["name"]}"')) + elif "TLC" in line: + fout.write(line.replace("TLC",f'"{self.configuration["system"]["structure1"]["tlc"]}"')) + else: + fout.write(line) + fin.close() + fout.close() def write_state( self, From 380d22e37eee894a0f41709a766278a7d5d97eda Mon Sep 17 00:00:00 2001 From: JohannesKarwou <72743318+JohannesKarwou@users.noreply.github.com> Date: Tue, 30 Aug 2022 22:19:11 +0200 Subject: [PATCH 4/9] small change in test --- transformato/tests/test_absolute.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/transformato/tests/test_absolute.py b/transformato/tests/test_absolute.py index 5ebb437e..ce07ba13 100644 --- a/transformato/tests/test_absolute.py +++ b/transformato/tests/test_absolute.py @@ -134,7 +134,7 @@ def test_perform_enstate_correction_asfe_system(): s1, mutation_list = create_asfe_system(configuration) - i = IntermediateStateFactory(system=s1, consecutive_runs=3, configuration=configuration) + i = IntermediateStateFactory(system=s1, configuration=configuration) perform_mutations( configuration=configuration, From c8b692a18b13a7e880737f0428b8a20c4587a0a5 Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Thu, 1 Sep 2022 12:55:53 +0200 Subject: [PATCH 5/9] pool all traj of intst1 if morge are available --- transformato/bin/perform_correction.py | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index 705ba576..7425612e 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -1,5 +1,6 @@ # general imports import os +import glob from endstate_correction.system import create_charmm_system, read_box from openmm.app import ( CharmmParameterSet, @@ -79,17 +80,13 @@ def gen_box(psf, crd): ######################################################## # ------------------- load samples ---------------------# - mm_samples = [] - if os.path.isfile(f"../{system_name}/intst1/run_1/lig_in_{env}.dcd"): - path_mm = f"../{system_name}/intst1/run_1/lig_in_{env}.dcd" - else: - path_mm = f"../{system_name}/intst1/lig_in_{env}.dcd" - - traj = mdtraj.load_dcd(path_mm, top=psf_file) - + files = glob.glob(f"../{system_name}/intst1/**/lig_in_{env}.dcd",recursive=True) + traj = mdtraj.load(files,top = psf_file) + if env == "waterbox": traj.image_molecules() + mm_samples = [] mm_samples.extend(traj.xyz * unit.nanometer) # NOTE: this is in nanometer! #################################################### From aea198880116ade7bfdde24c051da7e37fbc1d5c Mon Sep 17 00:00:00 2001 From: JohannesKarwou <72743318+JohannesKarwou@users.noreply.github.com> Date: Thu, 1 Sep 2022 16:35:46 +0200 Subject: [PATCH 6/9] Update perform_correction.py --- transformato/bin/perform_correction.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index 7425612e..600730dc 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -53,7 +53,7 @@ def gen_box(psf, crd): # load forcefiled files (ligand.str and toppar files) parms = () - file = f"../{system_name}/intst1/{tcl}" + file = f"../{system_name}/intst1/{tcl.lower()}" if os.path.isfile(f"{file}.str"): parms += (f"{file.lower()}.str",) else: From 6ea8e02c122abb6942d94d6bd0e2c6576b36088e Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Mon, 5 Sep 2022 18:31:17 +0200 Subject: [PATCH 7/9] adapting script to newest endstate_correction version --- transformato/bin/perform_correction.py | 40 ++++++-------------------- 1 file changed, 9 insertions(+), 31 deletions(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index 600730dc..b25b4a16 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -1,7 +1,7 @@ # general imports import os import glob -from endstate_correction.system import create_charmm_system, read_box +from endstate_correction.system import create_charmm_system, gen_box from openmm.app import ( CharmmParameterSet, CharmmPsfFile, @@ -14,33 +14,9 @@ from openmm import unit - # Variables will be generated by transformato system_name = NAMEofSYSTEM -tcl = TLC - - -def gen_box(psf, crd): - coords = crd.positions - - min_crds = [coords[0][0], coords[0][1], coords[0][2]] - max_crds = [coords[0][0], coords[0][1], coords[0][2]] - - for coord in coords: - min_crds[0] = min(min_crds[0], coord[0]) - min_crds[1] = min(min_crds[1], coord[1]) - min_crds[2] = min(min_crds[2], coord[2]) - max_crds[0] = max(max_crds[0], coord[0]) - max_crds[1] = max(max_crds[1], coord[1]) - max_crds[2] = max(max_crds[2], coord[2]) - - boxlx = max_crds[0] - min_crds[0] - boxly = max_crds[1] - min_crds[1] - boxlz = max_crds[2] - min_crds[2] - - psf.setBox(boxlx, boxly, boxlz) - return psf - +tlc = TLC for env in ["waterbox", "vacuum"]: @@ -53,7 +29,7 @@ def gen_box(psf, crd): # load forcefiled files (ligand.str and toppar files) parms = () - file = f"../{system_name}/intst1/{tcl.lower()}" + file = f"../{system_name}/intst1/{tlc.lower()}" if os.path.isfile(f"{file}.str"): parms += (f"{file.lower()}.str",) else: @@ -80,9 +56,9 @@ def gen_box(psf, crd): ######################################################## # ------------------- load samples ---------------------# - files = glob.glob(f"../{system_name}/intst1/**/lig_in_{env}.dcd",recursive=True) - traj = mdtraj.load(files,top = psf_file) - + files = glob.glob(f"../{system_name}/intst1/**/lig_in_{env}.dcd", recursive=True) + traj = mdtraj.load(files, top=psf_file) + if env == "waterbox": traj.image_molecules() @@ -103,4 +79,6 @@ def gen_box(psf, crd): ) r = perform_endstate_correction(fep_protocoll) - plot_endstate_correction_results(system_name, r, f"results_neq_unidirectional_{env}.png") + plot_endstate_correction_results( + system_name, r, f"results_neq_unidirectional_{env}.png" + ) From 61b8b229d0ef34a0dd2b5f90b3e22d6cfed76152 Mon Sep 17 00:00:00 2001 From: JohannesKarwou Date: Wed, 7 Sep 2022 15:51:20 +0200 Subject: [PATCH 8/9] change name of submi script --- transformato/bin/slurm_switching.sh | 16 ++++++++++++++++ transformato/bin/submit_switching_slurm.sh | 17 ----------------- transformato/state.py | 4 ++-- 3 files changed, 18 insertions(+), 19 deletions(-) create mode 100644 transformato/bin/slurm_switching.sh delete mode 100644 transformato/bin/submit_switching_slurm.sh diff --git a/transformato/bin/slurm_switching.sh b/transformato/bin/slurm_switching.sh new file mode 100644 index 00000000..e1be10b2 --- /dev/null +++ b/transformato/bin/slurm_switching.sh @@ -0,0 +1,16 @@ + + + + + + + + + + +export OPENMM_PRECISION='mixed' + +pwd; hostname; date + + +python perform_correction.py diff --git a/transformato/bin/submit_switching_slurm.sh b/transformato/bin/submit_switching_slurm.sh deleted file mode 100644 index 53fb4a94..00000000 --- a/transformato/bin/submit_switching_slurm.sh +++ /dev/null @@ -1,17 +0,0 @@ -#!/bin/bash -#SBATCH -p lgpu -#SBATCH --gres=gpu -#SBATCH --job-name=switching # Job name -#SBATCH --ntasks=1 # Run on a single CPU -#SBATCH --output=switching_%j.log # Standard output and error log - - -source /home/johannes/miniconda3/etc/profile.d/conda.sh -conda activate endstate - -export OPENMM_PRECISION='mixed' - -pwd; hostname; date - - -python perform_correction.py diff --git a/transformato/state.py b/transformato/state.py index 579bb145..880e829c 100644 --- a/transformato/state.py +++ b/transformato/state.py @@ -53,8 +53,8 @@ def endstate_correction(self): logger.info(f"Folder for endstate correction exist already") # copy submit script to the newly created folder - submit_switching_script_source = f"{self.configuration['bin_dir']}/submit_switching_slurm.sh" - submit_switchting_script_target = f"{self.path}/../endstate_correction/submit_switching_slurm.sh" + submit_switching_script_source = f"{self.configuration['bin_dir']}/slurm_switching.sh" + submit_switchting_script_target = f"{self.path}/../endstate_correction/slurm_switching.sh" shutil.copyfile(submit_switching_script_source, submit_switchting_script_target) # modify the perform_correction file from transformato bin and save it in the endstate_correcition folder From 9d47aa64a09a68fadab14691eb5f5241d538d1b7 Mon Sep 17 00:00:00 2001 From: JohannesKarwou <72743318+JohannesKarwou@users.noreply.github.com> Date: Fri, 9 Sep 2022 09:23:07 +0200 Subject: [PATCH 9/9] change paths in perform_correction --- transformato/bin/perform_correction.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/transformato/bin/perform_correction.py b/transformato/bin/perform_correction.py index b25b4a16..1888b0ca 100644 --- a/transformato/bin/perform_correction.py +++ b/transformato/bin/perform_correction.py @@ -20,8 +20,8 @@ for env in ["waterbox", "vacuum"]: - parameter_base = f"/site/raid3/johannes/free_solv_test/data" - # load the charmm specific files (psf, pdb, rtf, crd files) + + # load the charmm specific files (psf, rtf, crd files) psf_file = f"../{system_name}/intst1/lig_in_{env}.psf" psf = CharmmPsfFile(psf_file) pdb = PDBFile(f"../{system_name}/intst1/lig_in_{env}.pdb") @@ -31,10 +31,10 @@ parms = () file = f"../{system_name}/intst1/{tlc.lower()}" if os.path.isfile(f"{file}.str"): - parms += (f"{file.lower()}.str",) + parms += (f"{file}.str",) else: - parms += (f"{file.lower()}_g.rtf",) - parms += (f"{file.lower()}.prm",) + parms += (f"{file}_g.rtf",) + parms += (f"{file}.prm",) parms += (f"../toppar/top_all36_cgenff.rtf",) parms += (f"../toppar/par_all36_cgenff.prm",)