Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2013-2018 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 "ReferenceConfiguration.h"
23 : #include "ReferenceArguments.h"
24 : #include "ReferenceAtoms.h"
25 : #include "core/Value.h"
26 : #include "tools/OFile.h"
27 : #include "tools/PDB.h"
28 :
29 : namespace PLMD {
30 :
31 2135 : ReferenceConfigurationOptions::ReferenceConfigurationOptions( const std::string& type ):
32 2135 : tt(type)
33 : {
34 2135 : }
35 :
36 918 : bool ReferenceConfigurationOptions::usingFastOption() const {
37 918 : return (tt.find("-FAST")!=std::string::npos);
38 : }
39 :
40 7 : std::string ReferenceConfigurationOptions::getMultiRMSDType() const {
41 7 : plumed_assert( tt.find("MULTI-")!=std::string::npos );
42 7 : std::size_t dot=tt.find_first_of("MULTI-");
43 7 : return tt.substr(dot+6);
44 : }
45 :
46 2135 : ReferenceConfiguration::ReferenceConfiguration( const ReferenceConfigurationOptions& ro ):
47 2135 : name(ro.tt)
48 : // arg_ders(0),
49 : // atom_ders(0)
50 : {
51 2135 : weight=0.0;
52 2135 : }
53 :
54 2135 : ReferenceConfiguration::~ReferenceConfiguration()
55 : {
56 2135 : }
57 :
58 406 : std::string ReferenceConfiguration::getName() const {
59 406 : return name;
60 : }
61 :
62 393 : void ReferenceConfiguration::set( const PDB& pdb ) {
63 393 : line=pdb.getRemark();
64 393 : std::string ignore;
65 393 : if( parse("TYPE",ignore,true) ) {
66 12 : if(ignore!=name) error("mismatch for name");
67 : }
68 393 : if( !parse("WEIGHT",weight,true) ) weight=1.0;
69 393 : read( pdb );
70 393 : }
71 :
72 : // void ReferenceConfiguration::setNumberOfArguments( const unsigned& n ){
73 : // arg_ders.resize(n); tmparg.resize(n);
74 : // }
75 :
76 : // void ReferenceConfiguration::setNumberOfAtoms( const unsigned& n ){
77 : // atom_ders.resize(n);
78 : // }
79 :
80 : // bool ReferenceConfiguration::getVirial( Tensor& virout ) const {
81 : // if(virialWasSet) virout=virial;
82 : // return virialWasSet;
83 : // }
84 :
85 18 : void ReferenceConfiguration::parseFlag( const std::string&key, bool&t ) {
86 18 : Tools::parseFlag(line,key,t);
87 18 : }
88 :
89 0 : void ReferenceConfiguration::error(const std::string& msg) {
90 0 : plumed_merror("error reading reference configuration of type " + name + " : " + msg );
91 : }
92 :
93 1400 : void ReferenceConfiguration::setNamesAndAtomNumbers( const std::vector<AtomNumber>& numbers, const std::vector<std::string>& arg ) {
94 1400 : ReferenceAtoms* atoms=dynamic_cast<ReferenceAtoms*>( this );
95 1400 : if(!atoms) {
96 847 : plumed_massert( numbers.size()==0, "expecting no atomic positions");
97 : //setNumberOfAtoms( 0 );
98 : } else {
99 553 : atoms->setAtomNumbers( numbers );
100 : // setNumberOfAtoms( numbers.size() );
101 : }
102 : // Copy the arguments to the reference
103 1400 : ReferenceArguments* args=dynamic_cast<ReferenceArguments*>( this );
104 1400 : if(!args) {
105 551 : plumed_massert( arg.size()==0, "expecting no arguments");
106 : // setNumberOfArguments(0);
107 : } else {
108 849 : args->setArgumentNames( arg );
109 : // setNumberOfArguments( arg.size() );
110 : }
111 1400 : }
112 :
113 1494 : void ReferenceConfiguration::setReferenceConfig( const std::vector<Vector>& pos, const std::vector<double>& arg, const std::vector<double>& metric ) {
114 : // plumed_dbg_assert( pos.size()==atom_ders.size() && arg.size()==arg_ders.size() );
115 : // Copy the atomic positions to the reference
116 1494 : ReferenceAtoms* atoms=dynamic_cast<ReferenceAtoms*>( this );
117 1494 : if(!atoms) {
118 947 : plumed_massert( pos.size()==0, "expecting no atomic positions");
119 : } else {
120 1094 : std::vector<double> align_in( pos.size(), 1.0 ), displace_in( pos.size(), 1.0 );
121 1094 : atoms->setReferenceAtoms( pos, align_in, displace_in );
122 : }
123 : // Copy the arguments to the reference
124 1494 : ReferenceArguments* args=dynamic_cast<ReferenceArguments*>( this );
125 1494 : if(!args) {
126 547 : plumed_massert( arg.size()==0 && metric.size()==0, "expecting no arguments");
127 : } else {
128 947 : args->setReferenceArguments( arg, metric );
129 : }
130 1494 : }
131 :
132 306 : void ReferenceConfiguration::checkRead() {
133 306 : if(!line.empty()) {
134 0 : std::string msg="cannot understand the following words from the input line : ";
135 0 : for(unsigned i=0; i<line.size(); i++) {
136 0 : if(i>0) msg = msg + ", ";
137 0 : msg = msg + line[i];
138 : }
139 0 : error(msg);
140 : }
141 306 : }
142 :
143 298 : bool ReferenceConfiguration::isDirection() const {
144 298 : return ( name=="DIRECTION" );
145 : }
146 :
147 112776 : double ReferenceConfiguration::calculate( const std::vector<Vector>& pos, const Pbc& pbc, const std::vector<Value*>& vals,
148 : ReferenceValuePack& myder, const bool& squared ) const {
149 : // clearDerivatives();
150 112776 : std::vector<double> tmparg( vals.size() );
151 114683 : for(unsigned i=0; i<vals.size(); ++i) tmparg[i]=vals[i]->get();
152 114534 : return calc( pos, pbc, vals, tmparg, myder, squared );
153 : }
154 :
155 : // void ReferenceConfiguration::copyDerivatives( const ReferenceConfiguration* ref ){
156 : // plumed_dbg_assert( ref->atom_ders.size()==atom_ders.size() && ref->arg_ders.size()==arg_ders.size() );
157 : // for(unsigned i=0;i<atom_ders.size();++i) atom_ders[i]=ref->atom_ders[i];
158 : // for(unsigned i=0;i<arg_ders.size();++i) arg_ders[i]=ref->arg_ders[i];
159 : // virialWasSet=ref->virialWasSet; virial=ref->virial;
160 : // }
161 :
162 0 : void ReferenceConfiguration::print( OFile& ofile, const double& time, const double& weight, const double& lunits, const double& old_norm ) {
163 0 : ofile.printf("REMARK TIME=%f LOG_WEIGHT=%f OLD_NORM=%f\n",time, weight, old_norm );
164 0 : print( ofile, "%f", lunits ); // HARD CODED FORMAT HERE AS THIS IS FOR CHECKPOINT FILE
165 0 : }
166 :
167 206 : void ReferenceConfiguration::print( OFile& ofile, const std::string& fmt, const double& lunits ) {
168 206 : ofile.printf("REMARK TYPE=%s\n",getName().c_str() );
169 206 : ReferenceArguments* args=dynamic_cast<ReferenceArguments*>(this);
170 206 : if(args) args->printArguments( ofile, fmt );
171 206 : ReferenceAtoms* atoms=dynamic_cast<ReferenceAtoms*>(this);
172 206 : if(atoms) atoms->printAtoms( ofile, lunits );
173 206 : ofile.printf("END\n");
174 206 : }
175 :
176 9900 : double distance( const Pbc& pbc, const std::vector<Value*> & vals, ReferenceConfiguration* ref1, ReferenceConfiguration* ref2, const bool& squared ) {
177 : unsigned nder;
178 9900 : if( ref1->getReferencePositions().size()>0 ) nder=ref1->getReferenceArguments().size() + 3*ref1->getReferencePositions().size() + 9;
179 9900 : else nder=ref1->getReferenceArguments().size();
180 :
181 19800 : MultiValue myvals( 1, nder ); ReferenceValuePack myder( ref1->getReferenceArguments().size(), ref1->getReferencePositions().size(), myvals );
182 9900 : double dist1=ref1->calc( ref2->getReferencePositions(), pbc, vals, ref2->getReferenceArguments(), myder, squared );
183 : #ifndef NDEBUG
184 : // Check that A - B = B - A
185 : double dist2=ref2->calc( ref1->getReferencePositions(), pbc, vals, ref1->getReferenceArguments(), myder, squared );
186 : plumed_dbg_assert( fabs(dist1-dist2)<epsilon );
187 : #endif
188 19800 : return dist1;
189 : }
190 :
191 2523 : }
|