def test_should_not_wrap_a_group_in_a_group(self): rest = [restraints.SelectableRestraint()] grps = [restraints.RestraintGroup(rest, 1)] with mock.patch("meld.system.restraints.RestraintGroup.__init__", spec=True) as group_init: restraints.SelectivelyActiveCollection(grps, 1) self.assertEqual(group_init.call_count, 0)
def test_num_active_below_zero_should_raise(self): rest = [restraints.SelectableRestraint()] with self.assertRaises(RuntimeError): restraints.RestraintGroup(rest, -1)
def test_should_raise_on_empy_restraint_list(self): with self.assertRaises(RuntimeError): restraints.RestraintGroup([], 0)
def test_should_not_accept_non_selectable_restraint(self): rest = [restraints.NonSelectableRestraint()] with self.assertRaises(RuntimeError): restraints.RestraintGroup(rest, 1)
def test_should_accept_selectable_restraint(self): rest = [restraints.SelectableRestraint()] grp = restraints.RestraintGroup(rest, 1) self.assertEqual(len(grp.restraints), 1)
def test_num_active_should_be_set(self): rest = [restraints.SelectableRestraint()] grp = restraints.RestraintGroup(rest, 1) self.assertEqual(grp.num_active, 1)
def test_num_active_above_n_rest_should_raise(self): rest = [restraints.SelectableRestraint()] with self.assertRaises(RuntimeError): restraints.RestraintGroup(rest, 2)
def get_secondary_structure_restraints( system: interfaces.ISystem, scaler: restraints.RestraintScaler, ramp: Optional[restraints.TimeRamp] = None, torsion_force_constant: Optional[u.Quantity] = None, distance_force_constant: Optional[u.Quantity] = None, quadratic_cut: Optional[u.Quantity] = None, first_residue: Optional[indexing.ResidueIndex] = None, min_secondary_match: int = 4, filename: Optional[str] = None, content: Optional[str] = None, file: Optional[TextIO] = None, ) -> List[restraints.RestraintGroup]: """ Get a list of secondary structure restraints. Args: system: the system scaler: the force scaler ramp: the ramp, default is ConstantRamp torsion_force_constant: force constant for torsions, default 2.5e-2 kJ/mol/deg^2 distance_force_constant: force constant for distances, default 2.5 kJ/mol/nm^2 quadratic_cut: switch from quadratic to linear beyond this distance, default 0.2 nm min_secondary_match: minimum number of elements to match in secondary structure first_residue: residue at which these secondary structure restraints start filename: filename to open content: contents to process file: object to read from Returns A list of :class:`RestraintGroups` .. note:: Specify exactly one of filename, contents, file. """ torsion_force_constant = (2.5e-2 * u.kilojoule_per_mole / u.degree**2 if torsion_force_constant is None else torsion_force_constant) distance_force_constant = (2500 * u.kilojoule_per_mole / u.nanometer**2 if distance_force_constant is None else distance_force_constant) quadratic_cut = 0.2 * u.nanometer if quadratic_cut is None else quadratic_cut if ramp is None: ramp = restraints.ConstantRamp() if first_residue is None: first_residue = indexing.ResidueIndex(0) else: assert isinstance(first_residue, indexing.ResidueIndex) if min_secondary_match > 5: raise RuntimeError( "Minimum number of elements to match in secondary structure " "must be less than or equal to 5.") min_secondary_match = int(min_secondary_match) contents = _get_secondary_sequence(filename, content, file) groups = [] helices = _extract_secondary_runs(contents, "H", 5, min_secondary_match, first_residue) for helix in helices: rests: List[restraints.SelectableRestraint] = [] for index in range(helix.start + 1, helix.end - 1): phi = restraints.TorsionRestraint( system, scaler, ramp, system.index.atom(index - 1, "C"), system.index.atom(index, "N"), system.index.atom(index, "CA"), system.index.atom(index, "C"), -62.5 * u.degree, 17.5 * u.degree, torsion_force_constant, ) psi = restraints.TorsionRestraint( system, scaler, ramp, system.index.atom(index, "N"), system.index.atom(index, "CA"), system.index.atom(index, "C"), system.index.atom(index + 1, "N"), -42.5 * u.degree, 17.5 * u.degree, torsion_force_constant, ) rests.append(phi) rests.append(psi) d1 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(helix.start, "CA"), system.index.atom(helix.start + 3, "CA"), 0 * u.nanometer, 0.485 * u.nanometer, 0.561 * u.nanometer, 0.561 * u.nanometer + quadratic_cut, distance_force_constant, ) d2 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(helix.start + 1, "CA"), system.index.atom(helix.start + 4, "CA"), 0 * u.nanometer, 0.485 * u.nanometer, 0.561 * u.nanometer, 0.561 * u.nanometer + quadratic_cut, distance_force_constant, ) d3 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(helix.start, "CA"), system.index.atom(helix.start + 4, "CA"), 0 * u.nanometer, 0.581 * u.nanometer, 0.684 * u.nanometer, 0.684 * u.nanometer + quadratic_cut, distance_force_constant, ) rests.append(d1) rests.append(d2) rests.append(d3) group = restraints.RestraintGroup(rests, len(rests)) groups.append(group) extended = _extract_secondary_runs(contents, "E", 5, min_secondary_match, first_residue) for ext in extended: rests = [] for index in range(ext.start + 1, ext.end - 1): phi = restraints.TorsionRestraint( system, scaler, ramp, system.index.atom(index - 1, "C"), system.index.atom(index, "N"), system.index.atom(index, "CA"), system.index.atom(index, "C"), -117.5 * u.degree, 27.5 * u.degree, torsion_force_constant, ) psi = restraints.TorsionRestraint( system, scaler, ramp, system.index.atom(index, "N"), system.index.atom(index, "CA"), system.index.atom(index, "C"), system.index.atom(index + 1, "N"), 145 * u.degree, 25.0 * u.degree, torsion_force_constant, ) rests.append(phi) rests.append(psi) d1 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(ext.start, "CA"), system.index.atom(ext.start + 3, "CA"), 0 * u.nanometer, 0.785 * u.nanometer, 1.063 * u.nanometer, 1.063 * u.nanometer + quadratic_cut, distance_force_constant, ) d2 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(ext.start + 1, "CA"), system.index.atom(ext.start + 4, "CA"), 0 * u.nanometer, 0.785 * u.nanometer, 1.063 * u.nanometer, 1.063 * u.nanometer + quadratic_cut, distance_force_constant, ) d3 = restraints.DistanceRestraint( system, scaler, ramp, system.index.atom(ext.start, "CA"), system.index.atom(ext.start + 4, "CA"), 0 * u.nanometer, 1.086 * u.nanometer, 1.394 * u.nanometer, 1.394 * u.nanometer + quadratic_cut, distance_force_constant, ) rests.append(d1) rests.append(d2) rests.append(d3) group = restraints.RestraintGroup(rests, len(rests)) groups.append(group) return groups