The problem with MD prolongation
Hello, everyone!
I wrote the working MD script with gmxapi. However sometimes we need the prolongation of molecular dynamics and I would like to include it to my code. As I've found out the template of the command is right, but it didn't work. I caught the message 'Aborted (core dumped)'. I don't understand what is the reason of this problem. Please, help me to understand to resolve it.
Here is my code and the question is about the last strings in it.
import gmxapi as gmx
import os
input_dir= ''
input_dir = os.path.abspath(input_dir)
file1 = os.path.join(input_dir, 'WT.pdb')
print(f"Input file: {file1}")
assert os.path.isabs(file1)
assert os.path.exists(file1)
file2 = os.path.join(input_dir, "WT_processed.pdb")
pdb2gmx = gmx.commandline_operation('gmx',
arguments=['pdb2gmx', '-water', 'spce', '-ignh', 'yes', '-ff', 'amber99sb', '-missing', 'yes'],
input_files={'-f': file1},
output_files={"-o": file2,
'-p': os.path.join(input_dir, 'topol.top'),
'-i': os.path.join(input_dir, 'porse.itp')
})
pdb2gmx.run()
file3 = os.path.join(input_dir, "WT_newbox.pdb")
editconf = gmx.commandline_operation('gmx',
arguments=['editconf', '-c', 'yes', '-d', '1.0', '-bt', 'cubic'],
input_files={'-f': file2},
output_files={'-o': file3,
})
editconf.run()
file4 = os.path.join(input_dir, "WT_solv.pdb")
solvate = gmx.commandline_operation('gmx',
arguments=['solvate'],
input_files={'-cp': file3},
output_files={'-o': file4,
'-p': pdb2gmx.output.file['-p'],
})
solvate.run()
grompp = gmx.commandline_operation('gmx',
arguments=['grompp', '-maxwarn', '2'],
input_files={
'-f': os.path.join(input_dir, 'ions.mdp'),
'-p': solvate.output.file['-p'],
'-c': solvate.output.file['-o'],
# '-po': mdout_mdp,
},
output_files={
'-o': os.path.join(input_dir, 'ions.tpr')
})
grompp.run()
#grompp_file = grompp.output.file['-o'].result()
genion = gmx.commandline_operation('gmx',
arguments=['genion', '-pname', 'NA', '-nname', 'CL', '-neutral', 'yes'],
input_files={'-s': grompp.output.file['-o'],
'-p': solvate.output.file['-p']},
output_files={
'-o': os.path.join(input_dir, 'WT_solv_ions.pdb')},
# '-p': 'topol.top'},
stdin= 'SOL\n \n')
genion.run()
grompp2 = gmx.commandline_operation('gmx',
arguments=['grompp', '-maxwarn', '2'],
input_files={
'-f': os.path.join(input_dir, 'em.mdp'),
'-c': genion.output.file['-o'],
'-p': solvate.output.file['-p']
},
output_files={
'-o': os.path.join(input_dir, 'em.tpr'),
#'-p': os.path.join(input_dir, 'topol2.top')
})
grompp2.run()
err = grompp2.output.returncode.result()
simulation_input = gmx.read_tpr('em.tpr')
md = gmx.mdrun(input=simulation_input)
md.run()
md_path = md.output.directory.result()
energy = gmx.commandline_operation('gmx',
arguments=['energy'],
input_files={
'-f': os.path.join(md_path, 'ener.edr')},
output_files={
'-o': os.path.join(input_dir, 'em_potential.xvg')},
stdin='Potential\n \n')
energy.run()
grompp_nvt = gmx.commandline_operation('gmx',
arguments=['grompp', '-maxwarn', '2'],
input_files={
'-f': os.path.join(input_dir, 'nvt.mdp'),
'-c': os.path.join(md_path, 'confout.gro'),
'-r': os.path.join(md_path, 'confout.gro'),
'-p': os.path.join(input_dir, 'topol.top')
},
output_files={
'-o': os.path.join(input_dir, 'nvt.tpr'),
#'-p': os.path.join(input_dir, 'topol.top')
})
grompp_nvt.run()
simulation_input = gmx.read_tpr('nvt.tpr')
md = gmx.mdrun(input=simulation_input)
md.run()
md_path = md.output.directory.result()
energy_temp = gmx.commandline_operation('gmx',
arguments=['energy'],
input_files={
'-f': os.path.join(md_path, 'ener.edr')},
output_files={
'-o': os.path.join(input_dir, 'nvt_temp.xvg')},
stdin='Temperature\n \n')
energy_temp.run()
print(md.output.checkpoint.result())
grompp_npt = gmx.commandline_operation('gmx',
arguments=['grompp', '-maxwarn', '2'],
input_files={
'-f': os.path.join(input_dir, 'npt.mdp'),
'-c': os.path.join(md_path, 'confout.gro'),
'-r': os.path.join(md_path, 'confout.gro'),
'-t': md.output.checkpoint.result(),
'-p': os.path.join(input_dir, 'topol.top')
},
output_files={
'-o': os.path.join(input_dir, 'npt.tpr'),
#'-p': os.path.join(input_dir, 'topol.top')
})
grompp_npt.run()
simulation_input = gmx.read_tpr('npt.tpr')
md = gmx.mdrun(input=simulation_input)
md.run()
md_path = md.output.directory.result()
energy_press = gmx.commandline_operation('gmx',
arguments=['energy'],
input_files={
'-f': os.path.join(md_path, 'ener.edr')},
output_files={
'-o': os.path.join(input_dir, 'npt_press_dens.xvg')},
stdin='Pressure\nDensity\n \n')
energy_press.run()
grompp_md = gmx.commandline_operation('gmx',
arguments=['grompp', '-maxwarn', '2'],
input_files={
'-f': os.path.join(input_dir, 'md_5ns.mdp'),
'-c': os.path.join(md_path, 'confout.gro'),
'-r': os.path.join(md_path, 'confout.gro'),
'-t': md.output.checkpoint.result(),
'-p': os.path.join(input_dir, 'topol.top')
},
output_files={
'-o': os.path.join(input_dir, 'md_WT.tpr'),
#'-p': os.path.join(input_dir, 'topol.top')
})
grompp_md.run()
simulation_input = gmx.read_tpr('md_WT.tpr')
md = gmx.mdrun(input=simulation_input)
md.run()
md_path = md.output.directory.result()
trjconv = gmx.commandline_operation('gmx',
arguments=['trjconv', '-pbc', 'mol', '-center'],
input_files={
'-s': os.path.join(input_dir, 'md_WT.tpr'),
'-f': os.path.join(md_path, 'traj_comp.xtc')},
output_files={
'-o': os.path.join(input_dir, 'md_WT_protPBC.pdb')},
stdin='Protein\nProtein\n \n')
#trjconv.run()
convert_tpr = gmx.commandline_operation('gmx',
arguments=['convert-tpr', '-extend', '1000', ],
input_files={
'-s': os.path.join(input_dir, 'md_WT.tpr')
},
output_files={
'-o': os.path.join(input_dir, 'md_WT_cont_1ns.tpr')
})
convert_tpr.run()
err = convert_tpr.output.returncode.result()
if err != 0:
print(f"ERROR {err}: {convert_tpr.output.stderr.result()}")
else:
outfile = convert_tpr.output.file["-o"].result()
print(f"Output file path: {outfile}")
assert os.path.exists(outfile)
# Tutorials on gmxqpi and examples
# Here is about prolongation: https://gitlab.com/gromacs/gromacs/-/blob/main/python_packaging/gmxapi/test/test_mdrun.py
simulation_input = gmx.read_tpr('md_WT_cont_1ns.tpr')
md2 = gmx.mdrun(input=simulation_input, runtime_args={"-cpi": md.output.checkpoint})
md2.run()
Edited by Mariia Kotova