-
Notifications
You must be signed in to change notification settings - Fork 3
Expand file tree
/
Copy path_config.py
More file actions
2736 lines (2280 loc) · 100 KB
/
_config.py
File metadata and controls
2736 lines (2280 loc) · 100 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
######################################################################
# SOMD2: GPU accelerated alchemical free-energy engine.
#
# Copyright: 2023-2026
#
# Authors: The OpenBioSim Team <team@openbiosim.org>
#
# SOMD2 is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# SOMD2 is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with SOMD2. If not, see <http://www.gnu.org/licenses/>.
#####################################################################
"""
Configuration class for SOMD2 runner.
"""
__all__ = ["Config"]
from typing import Iterable as _Iterable
from openmm import Platform as _Platform
from pathlib import Path as _Path
import sire as _sr
from somd2 import _logger
# List of supported Sire platforms.
_sire_platforms = _sr.options.Platform.options()
# List of registered OpenMM platforms.
_omm_platforms = [
_Platform.getPlatform(x).getName().lower()
for x in range(0, _Platform.getNumPlatforms())
]
# List of available and supported platforms.
_platforms = ["auto"] + [x for x in _sire_platforms if x in _omm_platforms]
class Config:
"""
Class for storing a SOMD2 simulation configuration.
"""
# A dictionary of choices for options that support them. Here we inspect
# the Sire options module to get the valid choices. This allows us to be
# forwards compatible with new options.
_choices = {
"constraint": _sr.options.Constraint.options(),
"perturbable_constraint": _sr.options.PerturbableConstraint.options(),
"integrator": [
x
for x in _sr.options.Integrator.options()
if x not in ["auto", "verlet", "leapfrog"]
],
"cutoff_type": _sr.options.Cutoff.options(),
"platform": _platforms,
"lambda_schedule": [
"standard_morph",
"charge_scaled_morph",
"ring_break_morph",
"reverse_ring_break_morph",
],
"log_level": [level.lower() for level in _logger._core.levels],
"softcore_form": ["zacharias", "taylor"],
}
# A dictionary of nargs for the various options.
_nargs = {
"lambda_values": "+",
"lambda_energy": "+",
"rest2_scale": "+",
}
def __init__(
self,
log_level="info",
log_file="log.txt",
runtime="1 ns",
timestep="4 fs",
temperature="300 K",
pressure="1 atm",
surface_tension=None,
barostat_frequency=25,
integrator="langevin_middle",
cutoff_type="pme",
cutoff="7.5 A",
h_mass_factor=1.5,
hmr=True,
num_lambda=11,
lambda_values=None,
lambda_energy=None,
lambda_schedule="standard_morph",
charge_scale_factor=0.2,
swap_end_states=False,
coulomb_power=0.0,
shift_coulomb="1 A",
shift_delta="1.5 A",
restraints=None,
constraint="h_bonds",
perturbable_constraint="h_bonds_not_heavy_perturbed",
include_constrained_energies=False,
dynamic_constraints=True,
ghost_modifications=True,
fix_perturbable_zero_sigmas=True,
charge_difference=None,
coalchemical_restraint_dist=None,
com_reset_frequency=10,
minimise=True,
minimisation_constraints=False,
minimisation_errors=False,
equilibration_time="0 ps",
equilibration_timestep="2 fs",
equilibration_constraints=True,
energy_frequency="1 ps",
save_trajectories=True,
frame_frequency="100 ps",
save_velocities=False,
checkpoint_frequency="100 ps",
num_checkpoint_workers=None,
num_energy_neighbours=None,
null_energy="1e6 kcal/mol",
platform="auto",
max_threads=None,
max_gpus=None,
max_sire_threads=None,
opencl_platform_index=0,
oversubscription_factor=1,
replica_exchange=False,
randomise_velocities=False,
perturbed_system=None,
terminal_flip_frequency=None,
terminal_flip_angle=None,
terminal_flip_max_mobile_atoms=None,
gcmc=False,
gcmc_frequency=None,
gcmc_selection=None,
gcmc_excess_chemical_potential="-6.09 kcal/mol",
gcmc_standard_volume="30.543 A^3",
gcmc_num_waters=20,
gcmc_radius="4 A",
gcmc_bulk_sampling_probability=0.1,
gcmc_tolerance=0.0,
rest2_scale=1.0,
rest2_selection=None,
softcore_form="zacharias",
taylor_power=1,
output_directory="output",
restart=False,
use_backup=False,
write_config=True,
overwrite=False,
somd1_compatibility=False,
pert_file=None,
save_crash_report=False,
save_energy_components=False,
page_size=None,
timeout="300 s",
):
"""
Constructor.
Parameters
----------
runtime: str
Simulation length for each lambda window.
timestep: str
Integration time step.
temperature: str
Simulation temperature.
pressure: str
Simulation pressure. (Simulations will run in the NVT ensemble unless
a pressure is specified.)
surface_tension: str
Surface tension to use for NPT simulations with a membrane barostat.
barostat_frequency: int
The number of integration steps between barostat updates.
integrator: str
Integrator to use for simulation.
cutoff_type: str
Cutoff type to use for simulation.
cutoff: str
Non-bonded cutoff distance. Use "infinite" for no cutoff.
h_mass_factor: float
Factor by which to scale hydrogen masses.
hmr: bool
Whether to use hydrogen mass repartitioning. If False, then the masses
of the input system will be used. This can be useful if you have
already repartitioned the masses, or use a different repartitioning
scheme.
num_lambda: int
Number of lambda windows to use.
lambda_values: [float]
A list of lambda values. When specified, this takes precedence over
the 'num_lambda' option.
lambda_energy: [float]
A list of lambda values at which to output energy data. If not set,
then this will be set to the same as 'lambda_values', or the values
defined by 'num_lambda' if 'lambda_values' is not set.
lambda_schedule: str
Lambda schedule to use for alchemical free energy simulations.
charge_scale_factor: float
Factor by which to scale charges for charge scaled morph.
swap_end_states: bool
Whether to swap the end states of the alchemical system.
coulomb_power : float
Power to use for the soft-core Coulomb interaction. This is used
to soften the electrostatic interaction.
shift_coulomb : str
The soft-core shift-coulomb parameter. This is used to soften the
Coulomb interaction.
shift_delta : str
The soft-core shift-delta parameter. This is used to soften the
Lennard-Jones interaction.
restraints: sire.mm._MM.Restraints
A single set of restraints, or a list of sets of restraints that
will be applied to the atoms during the simulation.
constraint: str
Constraint type to use for non-perturbable molecules.
perturbable_constraint: str
Constraint type to use for perturbable molecules. If None, then
this will be set according to what is chosen for the non-perturbable
constraint.
include_constrained_energies: bool
Whether to include constrained energies in the potential.
dynamic_constraints: bool
Whether or not to update the length of constraints of perturbable
bonds with lambda. This defaults to True, meaning that changing
lambda will change any constraint on a perturbable bond to equal
to the value of r0 at that lambda value. If this is False, then
the constraint is set based on the current length.
ghost_modifications: bool
Whether to modify bonded terms between ghost atoms and the physical
system to avoid spurious coupling between the two, which can lead to
sampling of non-physical conformations. We implement the recommended
modifcations from https://pubs.acs.org/doi/10.1021/acs.jctc.0c01328
fix_perturbable_zero_sigmas: bool
Whether to prevent LJ sigma values being perturbed to zero.
charge_difference: int
The charge difference between the two end states. (Perturbed minus
reference.) If None, then alchemical ions will automatically be
added to keep the charge constant throughout the perturbation. If
specified, then the user defined value will take precedence. Note
the reference used for the charge difference is the same, regardless
of whether swap-end-states is set, i.e. the states are swapped after
the charge difference is calculated and alchemical ions are added.
coalchemical_restraint_dist: str
The minimum distance at which co-alchemical ions will be kept relative
to the centre of mass of the perturbable molecule in the system. This is
used to keep the co-alchemical ion in the bulk, preventing it from interacting
with the protein or ligand. If None, then no restraint will be applied.
Only functions for charge change perturbations.
com_reset_frequency: int
Frequency at which to reset the centre of mass of the system.
minimise: bool
Whether to minimise the system before simulation.
minimisation_constraints: bool
Whether to use constraints during minimisation. If False, then no
constraints will be used. If True, then the use of constraints will be
determined based on the value of 'equilibration_constraints'.
minimisation_errors: bool
Whether to raise an exception if a minimisation fails to converge.
equilibration_time: str
Time interval for equilibration. Only simulations starting from
scratch will be equilibrated.
equilibration_timestep: str
Equilibration timestep. (Can be different to simulation timestep.)
equilibration_constraints: bool
Whether to use constraints during equilibration. If False, then no constraints
will be used. If True, then the values specified by 'constraint' and
'perturbable_constraint' will be used.
energy_frequency: str
Frequency at which to output energy data. If running using 'replica_exchange',
then this will also be the frequency at which replica swaps are attempted.
When performing Grand Canonical Monte Carlo (GCMC) water insertions/deletions
via 'gcmc=True', this will also be the frequency at which GCMC moves are
attempted unless 'gcmc_frequency' is set.
save_trajectories: bool
Whether to save trajectory files
frame_frequency: str
Frequency at which to output trajectory frames.
save_velocities: bool
Whether to save velocities in trajectory frames.
checkpoint_frequency: str
Frequency at which to save checkpoint files, should be larger than
min(energy_frequency, frame_frequency). If zero, then no checkpointing
will be performed.
num_checkpoint_workers: int
The number of parallel workers to use when checkpointing during a replica
exchange simulation. By default, this is set to the number of concurrent
GPU contexts, i.e. the number of GPUs multiplied by the oversubscription
factor. The option can be used to reduce the number of workers, which
can be useful when the system size is large, i.e. when many large
trajectory files could be written simultaneously.
platform: str
Platform to run simulation on.
max_threads: int
Maximum number of CPU threads to use for simulation. (Default None, uses all available)
Does nothing if platform is set to CUDA.
max_gpus: int
Maximum number of GPUs to use for simulation (Default None, uses all available.)
Does nothing if platform is set to CPU.
max_sire_threads: int
Maximum number of CPU threads to use within Sire (e.g. for I/O operations).
(Default None, divides the total available threads between the number of
GPUs multiplied by the oversubscription factor.)
opencl_platform_index: int
The OpenCL platform index to use when multiple OpenCL implementations are
available on the system.
oversubscription_factor: int
The number of OpenMM contexts that can be run on a single GPU at the same time.
replica_exchange: bool
Whether to run replica exchange simulation. Currently this can only be used when
GPU resources are available.
randomise_velocities: bool
Whether to randomise velocities at the start of each replica exchange cycle
or following a terminal flip Monte Carlo move.
perturbed_system: str
The path to a stream file containing a Sire system for the equilibrated perturbed
end state (lambda = 1). This will be used as the starting conformation all lambda
windows > 0.5 when performing a replica exchange simulation.
terminal_flip_frequency: str
Frequency at which to attempt terminal ring flip Monte Carlo moves. If None
(the default), no terminal flip moves will be performed. When set, terminal
ring groups in perturbable molecules are detected automatically using Sire's
native connectivity. This must be a multiple of 'energy_frequency'.
terminal_flip_angle: str
Override the flip angle used for all terminal ring groups, e.g.
``"180 degrees"``. If None (the default), the angle is determined
automatically for each group from its geometry.
terminal_flip_max_mobile_atoms: int or None
Maximum number of mobile atoms allowed in a terminal ring group.
Groups with more mobile atoms than this threshold are skipped during
detection. Defaults to None (no limit).
gcmc: bool
Whether to perform Grand Canonical Monte Carlo (GCMC) water insertions/deletions.
gcmc_frequency: str
Frequency at which to attempt GCMC moves. If None, then this will be set to the
same as 'energy_frequency'. This must be a multiple of 'energy_frequency'.
gcmc_selection: str
A sire sslection string specifying the atoms that define the centre of geometry
of the GCMC sphere. If None, then GCMC moves will be attempted within the entire
simulation volume.
gcmc_excess_chemical_potential: str
The excess chemical potential of water in kcal/mol. The default value is calibrated
for the TIP3P water model. This can be calculated from the free energy of decoupling
a single water molecule from bulk.
gcmc_standard_volume: str
The standard volume of a water molecule in A^3. The default value is calibrated
from NPT simulation of TIP3P water.
gcmc_num_waters: int
The additional number of ghost water molecules to add to the system. These are
used as placeholders for GCMC insertion moves.
gcmc_radius: str
The radius of the GCMC sphere.
gcmc_bulk_sampling_probability: float
The probability of performing bulk GCMC moves, i.e. within the entire simulation
box rather than the GCMC sphere. These can be used to maintain a constant bulk
density, i.e. acting as a barostat. (This option has no affect when
'gcmc_selection=None'.)
gcmc_tolerance: float
The tolerance for the GCMC acceptance probability, i.e. the minimum probability
of acceptance for a move. This can be used to exclude low probability candidates
that can cause instabilities or crashes for the MD engine.
rest2_scale: float, list(float)
The scaling factor for Replica Exchange with Solute Tempering (REST) simulations.
This is the factor by which the temperature of the solute is scaled with respect to
the rest of the system. This can either be a single scaling factor, or a list of
scale factors for each lambda window. When a single scaling factor is used, then
the scale factor will be interpolated between a value of 1.0 in the end states,
and the value of 'rest2_scale' in intermediate lambda = 0.5 state. When multiple
values are used, then the number should match the number of lambda windows at which
energies are sampled.
rest2_selection: str
A sire selection string for atoms to include in the REST2 region in
addition to any perturbable molecules. For example, "molidx 0 and residx 0,1,2"
would select atoms from the first three residues of the first molecule. If None,
then all atoms within perturbable molecules will be included in the REST2 region.
When atoms within a perturbable molecule are included in the selection, then only
those atoms will be considered as part of the REST2 region. This allows REST2 to
be applied to protein mutations.
softcore_form: str
The soft-core potential form to use for alchemical interactions. This can be
either "zacharias" or "taylor". The default is "zacharias".
taylor_power: int
The power to use for the alpha term in the Taylor soft-core LJ expression,
i.e. sig6 = sigma^6 / (alpha^m * sigma^6 + r^6). Must be between 0 and 4.
The default is 1. Only used when softcore_form is "taylor".
output_directory: str
Path to a directory to store output files.
restart: bool
Whether to restart from a previous simulation using files found in 'output-directory'.
use_backup: bool
Whether to use backup files when restarting a simulation. If True, then
files from the last but one checkpoint will be used, rather than the most
recent checkpoint files. This can be useful if the most recent checkpoint
files are corrupted, or incomplete, e.g. you are recovering from a crash.
write_config: bool
Whether to write the configuration options to a YAML file in the output directory.
log_level: str
Log level to use.
log_file: str
Name of log file, will be saved in output directory.
overwrite: bool
Whether to overwrite files in the output directory, if files are detected and
this is false, SOMD2 will exit without overwriting.
somd1_compatibility: bool
Whether to run using a SOMD1 compatible perturbation.
pert_file: str
The path to a SOMD1 perturbation file to apply to the reference system.
When set, this will automatically set 'somd1_compatibility' to True.
save_crash_report: bool
Whether to save a crash report if the simulation crashes.
save_energy_components: bool
Whether to save the energy contribution for each force when checkpointing.
This is useful when debugging crashes.
page_size: int
The page size for trajectory handling in megabytes. If None, then Sire
will automatically set the page size.
timeout: str
Timeout for the minimiser and file lock.
num_energy_neighbours: int
The number of neighbouring windows to use when computing the energy
trajectory for the a given simulation lambda value. This can be
used to compute energies over a subset of windows, hence reducing
the cost of computing the energy trajectory. A value of 'null_energy'
will be added to the energy trajectory for the windows that are
omitted. If None, then all windows will be used.
null_energy: str
The energy value to use for lambda windows that are not
being computed as part of the energy trajectory.
"""
# Setup logger before doing anything else
self.log_level = log_level
self.log_file = log_file
self.output_directory = output_directory
self.runtime = runtime
self.temperature = temperature
self.pressure = pressure
self.surface_tension = surface_tension
self.barostat_frequency = barostat_frequency
self.integrator = integrator
self.cutoff_type = cutoff_type
self.cutoff = cutoff
self.h_mass_factor = h_mass_factor
self.hmr = hmr
self.timestep = timestep
self.num_lambda = num_lambda
self.lambda_values = lambda_values
self.lambda_energy = lambda_energy
self.lambda_schedule = lambda_schedule
self.charge_scale_factor = charge_scale_factor
self.swap_end_states = swap_end_states
self.coulomb_power = coulomb_power
self.shift_coulomb = shift_coulomb
self.shift_delta = shift_delta
self.restraints = restraints
self.constraint = constraint
self.perturbable_constraint = perturbable_constraint
self.include_constrained_energies = include_constrained_energies
self.dynamic_constraints = dynamic_constraints
self.ghost_modifications = ghost_modifications
self.fix_perturbable_zero_sigmas = fix_perturbable_zero_sigmas
self.charge_difference = charge_difference
self.coalchemical_restraint_dist = coalchemical_restraint_dist
self.com_reset_frequency = com_reset_frequency
self.minimise = minimise
self.minimisation_constraints = minimisation_constraints
self.minimisation_errors = minimisation_errors
self.equilibration_time = equilibration_time
self.equilibration_timestep = equilibration_timestep
self.equilibration_constraints = equilibration_constraints
self.energy_frequency = energy_frequency
self.save_trajectories = save_trajectories
self.frame_frequency = frame_frequency
self.save_velocities = save_velocities
self.checkpoint_frequency = checkpoint_frequency
self.num_checkpoint_workers = num_checkpoint_workers
self.platform = platform
self.max_threads = max_threads
self.max_gpus = max_gpus
self.max_sire_threads = max_sire_threads
self.opencl_platform_index = opencl_platform_index
self.oversubscription_factor = oversubscription_factor
self.replica_exchange = replica_exchange
self.randomise_velocities = randomise_velocities
self.perturbed_system = perturbed_system
self.terminal_flip_frequency = terminal_flip_frequency
self.terminal_flip_angle = terminal_flip_angle
self.terminal_flip_max_mobile_atoms = terminal_flip_max_mobile_atoms
self.gcmc = gcmc
self.gcmc_frequency = gcmc_frequency
self.gcmc_selection = gcmc_selection
self.gcmc_excess_chemical_potential = gcmc_excess_chemical_potential
self.gcmc_standard_volume = gcmc_standard_volume
self.gcmc_num_waters = gcmc_num_waters
self.gcmc_radius = gcmc_radius
self.gcmc_bulk_sampling_probability = gcmc_bulk_sampling_probability
self.gcmc_tolerance = gcmc_tolerance
self.rest2_scale = rest2_scale
self.rest2_selection = rest2_selection
self.restart = restart
self.use_backup = use_backup
self.softcore_form = softcore_form
self.taylor_power = taylor_power
self.somd1_compatibility = somd1_compatibility
self.pert_file = pert_file
self.save_crash_report = save_crash_report
self.save_energy_components = save_energy_components
self.timeout = timeout
self.num_energy_neighbours = num_energy_neighbours
self.null_energy = null_energy
self.page_size = page_size
self.write_config = write_config
self.overwrite = overwrite
def __str__(self):
"""Return a string representation of this object."""
# Get a dictionary representation of the object.
d = self.as_dict()
# Initialise the string.
string = "Config("
for k, v in d.items():
if isinstance(v, str):
string += f"{k.replace('', '')}='{v}', "
else:
string += f"{k.replace('', '')}={v}, "
# Remove the trailing comma and space.
string = string[:-2]
# Close the string.
string += ")"
return string
def __repr__(self):
"""Return a string representation of this object."""
return self.__str__()
def __eq__(self, other):
"""Equality operator."""
return self.as_dict() == other.as_dict()
@staticmethod
def from_yaml(path):
"""
Create a Config object from a YAML file.
Parameters
----------
path: str
Path to YAML file.
"""
from ..io import yaml_to_dict as _yaml_to_dict
d = _yaml_to_dict(path)
return Config(**d)
def as_dict(self, sire_compatible=False):
"""Convert config object to dictionary
Parameters
----------
sire_compatible: bool
Whether to convert to a dictionary compatible with Sire,
this simply converts any options with a value of None to a
boolean with the value False.
"""
from pathlib import Path as _Path
from sire.cas import LambdaSchedule as _LambdaSchedule
d = {}
for attr, value in self.__dict__.items():
if attr.startswith("_extra") or attr.startswith("extra"):
continue
attr_l = attr[1:]
if isinstance(value, _Path):
d[attr_l] = str(value)
else:
try:
d[attr_l] = value.to_string()
except AttributeError:
d[attr_l] = value
if value is None and sire_compatible:
d[attr_l] = False
# Handle the lambda schedule separately so that we can use simplified
# keyword options.
# A keyword exists for this lambda schedule.
if self._lambda_schedule_name is not None:
d["lambda_schedule"] = self._lambda_schedule_name
# Try to match the lambda schedule to a known schedule, if not then convert to hex.
else:
if self.lambda_schedule == _LambdaSchedule.standard_morph():
d["lambda_schedule"] = "standard_morph"
elif self.lambda_schedule == _LambdaSchedule.charge_scaled_morph(
self._charge_scale_factor
):
d["lambda_schedule"] = "charge_scaled_morph"
else:
d["lambda_schedule"] = self._to_hex(self.lambda_schedule)
# Serialise restraints.
if self.restraints is not None:
d["restraints"] = [self._to_hex(restraint) for restraint in self.restraints]
# Use the path for the perturbed_system option, since the system
# isn't serializable.
if (
self.perturbed_system is not None
and self._perturbed_system_file is not None
):
d["perturbed_system"] = str(self._perturbed_system_file)
d.pop("perturbed_system_file", None)
return d
@property
def runtime(self):
return self._runtime
@runtime.setter
def runtime(self, runtime):
if not isinstance(runtime, str):
raise TypeError("'runtime' must be of type 'str'")
from sire.units import picosecond
try:
t = _sr.u(runtime)
except:
raise ValueError(
f"Unable to parse 'runtime' as a Sire GeneralUnit: {runtime}"
)
if t.value() != 0 and not t.has_same_units(picosecond):
raise ValueError("'runtime' units are invalid.")
if t.value() == 0:
_logger.warning(
"Runtime is zero - simulation will not run. Set 'runtime' to a non-zero value."
)
self._runtime = t
@property
def temperature(self):
return self._temperature
@temperature.setter
def temperature(self, temperature):
if not isinstance(temperature, str):
raise TypeError("'temperature' must be of type 'str'")
from sire.units import kelvin
try:
t = _sr.u(temperature)
except:
raise ValueError(
f"Unable to parse 'temperature' as a Sire GeneralUnit: {temperature}"
)
if not t.has_same_units(kelvin):
raise ValueError("'temperature' units are invalid.")
self._temperature = t
@property
def pressure(self):
return self._pressure
@pressure.setter
def pressure(self, pressure):
if pressure is not None and not isinstance(pressure, str):
raise TypeError("'pressure' must be of type 'str'")
from sire.units import atm
if pressure is not None:
try:
p = _sr.u(pressure)
except:
# Handle special case of pressure = "none"
if pressure.lower().replace(" ", "") == "none":
self._pressure = None
return
raise ValueError(
f"Unable to parse 'pressure' as a Sire GeneralUnit: {pressure}"
)
if not p.has_same_units(atm):
raise ValueError("'pressure' units are invalid.")
self._pressure = p
else:
self._pressure = pressure
@property
def barostat_frequency(self):
return self._barostat_frequency
@barostat_frequency.setter
def barostat_frequency(self, barostat_frequency):
if not isinstance(barostat_frequency, int):
raise TypeError("'barostat_frequency' must be of type 'int'")
if barostat_frequency <= 0:
raise ValueError("'barostat_frequency' must be a positive integer")
self._barostat_frequency = barostat_frequency
@property
def surface_tension(self):
return self._surface_tension
@surface_tension.setter
def surface_tension(self, surface_tension):
if surface_tension is not None and not isinstance(surface_tension, str):
raise TypeError("'surface_tension' must be of type 'str'")
from sire.units import atm, angstrom
if surface_tension is not None:
try:
st = _sr.u(surface_tension)
except:
raise ValueError(
f"Unable to parse 'surface_tension' as a Sire GeneralUnit: {surface_tension}"
)
# Make sure we can handle a value of zero.
if st == 0:
st = 0 * atm * angstrom
elif not st.has_same_units(atm * angstrom):
raise ValueError("'surface_tension' units are invalid.")
self._surface_tension = st
else:
self._surface_tension = surface_tension
@property
def integrator(self):
return self._integrator
@integrator.setter
def integrator(self, integrator):
if not isinstance(integrator, str):
raise TypeError("'integrator' must be of type 'str'")
integrator = integrator.lower().replace(" ", "")
if integrator not in self._choices["integrator"]:
raise ValueError(
f"Integrator not recognised. Valid integrators are: {', '.join(self._choices['integrator'])}"
)
self._integrator = integrator
@property
def cutoff_type(self):
return self._cutoff_type
@cutoff_type.setter
def cutoff_type(self, cutoff_type):
if not isinstance(cutoff_type, str):
raise TypeError("'cutoff_type' must be of type 'str'")
cutoff_type = cutoff_type.lower().replace(" ", "")
if cutoff_type not in self._choices["cutoff_type"]:
raise ValueError(
f"Cutoff type not recognised. Valid cutoff types are: {', '.join(self._choices['cutoff_type'])}"
)
self._cutoff_type = cutoff_type
@property
def cutoff(self):
return self._cutoff
@cutoff.setter
def cutoff(self, cutoff):
if not isinstance(cutoff, str):
raise TypeError("'cutoff' must be of type 'str'")
from sire.units import angstrom
if cutoff is not None:
# Handle special case of cutoff = "infinite"
if cutoff.lower().replace(" ", "") == "infinite":
self._cutoff = "infinite"
else:
try:
c = _sr.u(cutoff)
except:
raise ValueError(
f"Unable to parse 'cutoff' as a Sire GeneralUnit: {cutoff}"
)
if not c.has_same_units(angstrom):
raise ValueError("'cutoff' units are invalid.")
self._cutoff = c
else:
self._cutoff = cutoff
@property
def h_mass_factor(self):
return self._h_mass_factor
@h_mass_factor.setter
def h_mass_factor(self, h_mass_factor):
if not isinstance(h_mass_factor, float):
try:
h_mass_factor = float(h_mass_factor)
except Exception:
raise ValueError("'h_mass_factor' must be a float")
if h_mass_factor < 1.0:
_logger.warning(
"Requested hydrogen mass repartitioning factor is less than 1.0. "
"This will result in a reduction of the mass of hydrogen atoms, "
"and will likely lead to undesired simulation behaviour."
)
self._h_mass_factor = h_mass_factor
@property
def hmr(self):
return self._hmr
@hmr.setter
def hmr(self, hmr):
if not isinstance(hmr, bool):
raise ValueError("'hmr' must be of type 'bool'")
self._hmr = hmr
@property
def timestep(self):
return self._timestep
@timestep.setter
def timestep(self, timestep):
if not isinstance(timestep, str):
raise TypeError("'timestep' must be of type 'str'")
from sire.units import femtosecond
try:
t = _sr.u(timestep)
except:
raise ValueError(
f"Unable to parse 'timestep' as a Sire GeneralUnit: {timestep}"
)
if t.value() != 0 and not t.has_same_units(femtosecond):
raise ValueError("'timestep' units are invalid.")
if t.value() == 0:
_logger.warning(
"Timestep is zero - simulation will not run. Set 'timestep' to a non-zero value."
)
if t > _sr.u("2fs") and self.h_mass_factor <= 1.0:
_logger.warning("Timestep is large - consider repartitioning hydrogen mass")
self._timestep = t
@property
def num_lambda(self):
return self._num_lambda
@num_lambda.setter
def num_lambda(self, num_lambda):
if num_lambda is not None:
if not isinstance(num_lambda, int):
raise ValueError("'num_lambda' must be an integer")
self._num_lambda = num_lambda
@property
def lambda_values(self):
return self._lambda_values
@lambda_values.setter
def lambda_values(self, lambda_values):
if lambda_values is not None:
if not isinstance(lambda_values, _Iterable):
raise ValueError("'lambda_values' must be an iterable")
try:
lambda_values = [float(x) for x in lambda_values]
except:
raise ValueError("'lambda_values' must be an iterable of floats")
if not all(0 <= x <= 1 for x in lambda_values):
raise ValueError(
"All entries in 'lambda_values' must be between 0 and 1"
)
# Round to 5dp.
lambda_values = [round(x, 5) for x in lambda_values]
self._num_lambda = len(lambda_values)