Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-2023 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #include "MultiColvarBase.h"
23 : #include "AtomValuePack.h"
24 : #include "tools/Torsion.h"
25 : #include "core/ActionRegister.h"
26 :
27 : #include <string>
28 : #include <cmath>
29 :
30 : namespace PLMD {
31 : namespace multicolvar {
32 :
33 : //+PLUMEDOC MCOLVAR TORSIONS
34 : /*
35 : Calculate whether or not a set of torsional angles are within a particular range.
36 :
37 : \par Examples
38 :
39 : The following provides an example of the input for the TORSIONS command
40 :
41 : \plumedfile
42 : TORSIONS ...
43 : ATOMS1=168,170,172,188
44 : ATOMS2=170,172,188,190
45 : ATOMS3=188,190,192,230
46 : BETWEEN={GAUSSIAN LOWER=0 UPPER=pi SMEAR=0.1}
47 : LABEL=ab
48 : ... TORSIONS
49 : PRINT ARG=ab.* FILE=colvar STRIDE=10
50 : \endplumedfile
51 :
52 : Writing out the atoms involved in all the torsion angles in this way can be rather tedious. Thankfully if you are working with protein you
53 : can avoid this by using the \ref MOLINFO command. PLUMED uses the pdb file that you provide to this command to learn
54 : about the topology of the protein molecule. This means that you can specify torsion angles using the following syntax:
55 :
56 : \plumedfile
57 : #SETTINGS MOLFILE=regtest/basic/rt32/helix.pdb
58 : MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
59 : TORSIONS ...
60 : ATOMS1=@phi-3
61 : ATOMS2=@psi-3
62 : ATOMS3=@phi-4
63 : BETWEEN={GAUSSIAN LOWER=0 UPPER=pi SMEAR=0.1}
64 : LABEL=ab
65 : ... TORSIONS
66 : PRINT ARG=ab.* FILE=colvar STRIDE=10
67 : \endplumedfile
68 :
69 : Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
70 : Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the fourth residue of the protein.
71 :
72 :
73 : */
74 : //+ENDPLUMEDOC
75 :
76 : class Torsions : public MultiColvarBase {
77 : public:
78 : static void registerKeywords( Keywords& keys );
79 : explicit Torsions(const ActionOptions&);
80 : double compute( const unsigned& tindex, AtomValuePack& myatoms ) const override;
81 5 : bool isPeriodic() override { return true; }
82 2 : void retrieveDomain( std::string& min, std::string& max ) override { min="-pi"; max="pi"; }
83 : };
84 :
85 10423 : PLUMED_REGISTER_ACTION(Torsions,"TORSIONS")
86 :
87 3 : void Torsions::registerKeywords( Keywords& keys ) {
88 3 : MultiColvarBase::registerKeywords( keys );
89 6 : keys.add("numbered","ATOMS","the atoms involved in each of the torsion angles you wish to calculate. "
90 : "Keywords like ATOMS1, ATOMS2, ATOMS3,... should be listed and one torsion will be "
91 : "calculated for each ATOM keyword you specify (all ATOM keywords should "
92 : "provide the indices of four atoms). The eventual number of quantities calculated by this "
93 : "action will depend on what functions of the distribution you choose to calculate.");
94 6 : keys.reset_style("ATOMS","atoms");
95 6 : keys.use("BETWEEN"); keys.use("HISTOGRAM");
96 3 : }
97 :
98 2 : Torsions::Torsions(const ActionOptions&ao):
99 : Action(ao),
100 2 : MultiColvarBase(ao)
101 : {
102 : // Read in the atoms
103 2 : int natoms=4; std::vector<AtomNumber> all_atoms;
104 2 : readAtomsLikeKeyword( "ATOMS", natoms, all_atoms );
105 2 : setupMultiColvarBase( all_atoms );
106 2 : std::vector<bool> catom_ind(4, false);
107 2 : catom_ind[1]=catom_ind[2]=true;
108 2 : setAtomsForCentralAtom( catom_ind );
109 : // Read in the vessels
110 2 : readVesselKeywords();
111 : // And check everything has been read in correctly
112 2 : checkRead();
113 2 : }
114 :
115 36 : double Torsions::compute( const unsigned& tindex, AtomValuePack& myatoms ) const {
116 36 : Vector d0,d1,d2;
117 36 : d0=getSeparation(myatoms.getPosition(1),myatoms.getPosition(0));
118 36 : d1=getSeparation(myatoms.getPosition(2),myatoms.getPosition(1));
119 36 : d2=getSeparation(myatoms.getPosition(3),myatoms.getPosition(2));
120 :
121 36 : Vector dd0,dd1,dd2; PLMD::Torsion t;
122 36 : double value = t.compute(d0,d1,d2,dd0,dd1,dd2);
123 :
124 36 : addAtomDerivatives(1,0,dd0,myatoms);
125 36 : addAtomDerivatives(1,1,dd1-dd0,myatoms);
126 36 : addAtomDerivatives(1,2,dd2-dd1,myatoms);
127 36 : addAtomDerivatives(1,3,-dd2,myatoms);
128 :
129 36 : myatoms.addBoxDerivatives (1, -(extProduct(d0,dd0)+extProduct(d1,dd1)+extProduct(d2,dd2)));
130 :
131 36 : return value;
132 : }
133 :
134 : }
135 : }
|
