The general procedure to set up a quantum refinement job consists of the following steps:
Build a model.
phenix.ready_set add_h_to_water=True
can be useful for this purpose.Prepare restraint files for unknown residues and ligands. The script qref_prep.py
will tell you if there are any missing restraint files.
phenix.ready_set
and phenix.elbow
.Prepare syst1
files that define the QM regions.
syst1
. It there are multiple QM regions, the names should be syst11
, syst12
, etc.sort_atoms = False
in the input to Phenix
should
ensure that the ordering in the input model is preserved, we have
encountered instances where this is not adhered to. Thus the script sort_pdb.py
is also provided with QRef
; this script takes as input the name of a PDB file and outputs a new PDB file with the suffix _sorted.pdb
where the atoms are sorted in the same order as that Phenix
uses internally. You should use this file as the input model for refinement, as well as the reference when defining the syst1
files.syst1
files allows for multiple atoms or intervals of atoms to be specified on a single line, where ,
or blank
works as delimiters and-
is used to indicate an interval.#
and !
can be used to include comments in the syst1
files.syst1
files
will indicate that this is a link atom, i.e. it will be replaced by a
hydrogen at the appropriate position in the QM calculation.Run qref_prep.py <model>_sorted.pdb
to generate qref.dat
(a file containing settings for the QR interface), as well as PDB files describing the QM regions.
The junctfactor
file needs to be present in the same directory as where qref_prep.py
is run. The junctfactor
file contain ideal QM distances for the junctfactor
file is to be used this can be specified with the -j
or --junctfactor
option. Our standard file is available on GitHub.
The theory used for the ideal QM distance can be changed with the -l
or --ltype
option. Default is type 12. Current options are:
There needs to be a parametrisation in the junctfactor
file for the bond one intends to cleave; it is recommended that the user inspects the junctfactor
file
to verify that there is support to cleave the intended bond type. In
the case parametrisation is lacking another selection for the QM system
(and in particular where the link between QM and MM occurs) needs to be
made or appropriate parametrisation added to the junctfactor
file.
Parametrisation
is done by simply optimising a small model of the QM system with the QM
method and then insert the optimised bond length in the juncfactor
file:
CYM CB CA H1 # Cym CB junction
7 1.1137d0 # MeS-, Cs, BP/RI, 6-31G*, 29/10-02
9 1.1307d0 # MeS-, Cs, BP/RI, def2-SV(P), UR 24/6-08
10 1.1306d0 # MeS-, Cs, PBE/RI, def2-SV(P), H in plane, UR 22/6-11
11 1.1201d0 # MeS-, Cs, B3LYP/RI, def2-SV(P), UR 18/11-09
12 1.1238d0 # MeS-, Cs, TPSS/RI, def2-SV(P), UR 29/12-09
The first line contains the residue name (CYM), the QM atom bound to the
link atom (CB) and the link atom (CA), i.e. the atom that is converted
to a H atom in the QM calculation. The fourth entry is the desired atom
type of the link atom in the MM force filed. It is not used in this
implementation.
The coming lines give the theory (ltype) and the optiised bond length in
Å. Characters after "#" are comments and are ignored. The entry ends
with an empty line.
Ideally only the input model is needed as an argument for qref_prep.py
. If there was a need to prepare restraint files for novel residues or ligands in point 2 above, qref_prep.py
needs to be made aware of these. This can be achieved with the -c
or --cif
option.
The output from qref_prep.py
should be (1 + 2 nsyst1) files as well as selection strings:
qref.dat
,
which contains the settings for the QR interface. This file can be
changed manually and it is a good idea to inspect that the value for orca_binary
is the correct path for the actual Orca binary file (qref_prep.py
tries to locate this file automatically but may sometimes fail).qm_i_c.pdb
, which is the model used to calculate the MM1 energy.qm_i_h.pdb
, which is the model used to calculate the QM1 energy.The output PDB files can be used to inspect that the QM selection is proper.
phenix.refine
or phenix.real_space_refine
, see point 6 below.All available options for qref_prep.py
can be seen through the -h
or --help
option.
The next step is to prepare input files for Orca.
qm_i.inp
, where i
refers to the i:th QM region; counting starts at 1.Orca
then becomes
! TPSS D4 DEF2-SV(P)
.
! ENGRAD
qm_i_h.pdb
) use
*pdbfile <charge> <multiplicity> qm_i_h.pdb
i
refers to the i:th QM region.examples
folder on GitHub.Now, you can start the quantum refinement. However, it is recommended to first create an empty file named qm.lock
(this disables QR; touch qm.lock), then start a refinement so that an .eff
file is obtained:
phenix.refine <model>_sorted.pdb <model>.mtz
The .eff file contains all the Phenix
refinement settings for the current experimental data and model). Rename this file to phenix.params
and edit it to make sure the following options are set:
phenix.refine
):refinement.pdb_interpretation.restraints_library.cdl = False
refinement.pdb_interpretation.restraints_library.mcl = False
refinement.pdb_interpretation.restraints_library.cis_pro_eh99 = False
refinement.pdb_interpretation.sort_atoms = False
refinement.refine.strategy
= *individual_sites individual_sites_real_space rigid_body
*individual_adp group_adp tls occupancies group_anomalous
refinement.refine.sites.individual = <reciprocal selection string>
refinement.hydrogens.refine = *individual riding Auto
refinement.main.nqh_flips = False
phenix.real_space_refine
):refinement.run = *minimization_global rigid_body local_grid_search morphing simulated_annealing adp occupancy nqh_flips
pdb_interpretation.restraints_library.cdl = False
pdb_interpretation.restraints_library.mdl = False
pdb_interpretation.restraints_library.cis_pro_eh99 = False
pdb_interpretation.sort_atoms = False
pdb_interpretation.secondary_structure = False
pdb_interpretation.reference_coordinate_restraints.enabled = True
pdb_interpretation.reference_coordinate_restraints.selection = <real space selection string>
pdb_interpretation.reference_coordinate_restraints.sigma = 0.01
To run the quantum refinement job, delete the qm.lock
file and then execute either phenix.refine phenix.params
or phenix.real_space_refine phenix.params
.
If there is a need to restart the job with different settings for Phenix
, make sure to delete the file settings.pickle
.
phenix.refine phenix.params --overwrite
or
phenix_real_space.refine phenix.params --overwrite
If the job is run in parallel, you need to load modules for openMPI
If you need to do a restart, do
rm *.lock
! TPSS D4 DEF2-SV(P)
! ENGRAD
%PAL NPROCS 4 END
*pdbfile +1 1 qm_1_h.pdb
A somewhat more hands-on description of the installation than described on GitHub.