Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-2019 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 : using namespace std;
31 :
32 : namespace PLMD {
33 : namespace multicolvar {
34 :
35 : //+PLUMEDOC MCOLVAR TORSIONS
36 : /*
37 : Calculate whether or not a set of torsional angles are within a particular range.
38 :
39 : \par Examples
40 :
41 : The following provides an example of the input for the torsions command
42 :
43 : \plumedfile
44 : TORSIONS ...
45 : ATOMS1=168,170,172,188
46 : ATOMS2=170,172,188,190
47 : ATOMS3=188,190,192,230
48 : LABEL=ab
49 : ... TORSIONS
50 : PRINT ARG=ab.* FILE=colvar STRIDE=10
51 : \endplumedfile
52 :
53 : Writing out the atoms involved in all the torsions in this way can be rather tedious. Thankfully if you are working with protein you
54 : can avoid this by using the \ref MOLINFO command. PLUMED uses the pdb file that you provide to this command to learn
55 : about the topology of the protein molecule. This means that you can specify torsion angles using the following syntax:
56 :
57 : \plumedfile
58 : MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
59 : TORSIONS ...
60 : ATOMS1=@phi-3
61 : ATOMS2=@psi-3
62 : ATOMS3=@phi-4
63 : LABEL=ab
64 : ... TORSIONS
65 : PRINT ARG=ab FILE=colvar STRIDE=10
66 : \endplumedfile
67 :
68 : Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
69 : Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the 4th residue of the protein.
70 :
71 :
72 : */
73 : //+ENDPLUMEDOC
74 :
75 0 : class Torsions : public MultiColvarBase {
76 : public:
77 : static void registerKeywords( Keywords& keys );
78 : explicit Torsions(const ActionOptions&);
79 : virtual double compute( const unsigned& tindex, AtomValuePack& myatoms ) const ;
80 0 : bool isPeriodic() { return true; }
81 0 : void retrieveDomain( std::string& min, std::string& max ) { min="-pi"; max="pi"; }
82 : };
83 :
84 6452 : PLUMED_REGISTER_ACTION(Torsions,"TORSIONS")
85 :
86 1 : void Torsions::registerKeywords( Keywords& keys ) {
87 1 : MultiColvarBase::registerKeywords( keys );
88 4 : keys.reserve("numbered","ATOMS","the atoms involved in each of the torsion angles you wish to calculate. "
89 : "Keywords like ATOMS1, ATOMS2, ATOMS3,... should be listed and one torsion will be "
90 : "calculated for each ATOM keyword you specify (all ATOM keywords should "
91 : "provide the indices of four atoms). The eventual number of quantities calculated by this "
92 : "action will depend on what functions of the distribution you choose to calculate.");
93 3 : keys.reset_style("ATOMS","atoms");
94 3 : keys.use("BETWEEN"); keys.use("HISTOGRAM");
95 1 : }
96 :
97 0 : Torsions::Torsions(const ActionOptions&ao):
98 : Action(ao),
99 0 : MultiColvarBase(ao)
100 : {
101 : // Read in the atoms
102 0 : int natoms=4; std::vector<AtomNumber> all_atoms;
103 0 : readAtomsLikeKeyword( "ATOMS", natoms, all_atoms );
104 0 : setupMultiColvarBase( all_atoms );
105 0 : std::vector<bool> catom_ind(4, false);
106 : catom_ind[1]=catom_ind[2]=true;
107 0 : setAtomsForCentralAtom( catom_ind );
108 : // Read in the vessels
109 0 : readVesselKeywords();
110 : // And check everything has been read in correctly
111 0 : checkRead();
112 0 : }
113 :
114 0 : double Torsions::compute( const unsigned& tindex, AtomValuePack& myatoms ) const {
115 0 : Vector d0,d1,d2;
116 0 : d0=getSeparation(myatoms.getPosition(1),myatoms.getPosition(0));
117 0 : d1=getSeparation(myatoms.getPosition(2),myatoms.getPosition(1));
118 0 : d2=getSeparation(myatoms.getPosition(3),myatoms.getPosition(2));
119 :
120 0 : Vector dd0,dd1,dd2; PLMD::Torsion t;
121 0 : double value = t.compute(d0,d1,d2,dd0,dd1,dd2);
122 :
123 0 : addAtomDerivatives(1,0,dd0,myatoms);
124 0 : addAtomDerivatives(1,1,dd1-dd0,myatoms);
125 0 : addAtomDerivatives(1,2,dd2-dd1,myatoms);
126 0 : addAtomDerivatives(1,3,-dd2,myatoms);
127 :
128 0 : myatoms.addBoxDerivatives (1, -(extProduct(d0,dd0)+extProduct(d1,dd1)+extProduct(d2,dd2)));
129 :
130 0 : return value;
131 : }
132 :
133 : }
134 4839 : }
|