Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2017-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 "ActionToPutData.h"
23 : #include "ActionRegister.h"
24 : #include "PlumedMain.h"
25 : #include "ActionSet.h"
26 :
27 : //+PLUMEDOC ANALYSIS PUT
28 : /*
29 : Pass data into PLUMED
30 :
31 : The PUT command transfers the data from a void pointer passed to PLUMED from the calling code to a PLMD::Value. The calling
32 : object knows the shapes of the variables it is passing, so if you want to pass a 3x3 matrix from the MD code to PLUMED, you create the space to do so as follows:
33 :
34 : ```c++
35 : plumed.cmd("readInputLine n: PUT SHAPE=3,3 UNIT=length PERIODIC=NO");
36 : ```
37 :
38 : This command then creates a PLMD::Value called `n` that you can refer to later in your PLUMED input file. To transfer data from the void pointer called val into the PLMD::Value
39 : called `n`, you would then use the following command:
40 :
41 : ```c++
42 : plumed.cmd("setValue n", val);
43 : ```
44 :
45 : Notice also that if you expect PLUMED to try to apply forces on `n`, you can pass a void pointer called `force` to get the forces that PLUMED has applied on the elements of n as follows:
46 :
47 : ```c++
48 : plumed.cmd("setValueForces n", force);
49 : ```
50 :
51 : Within the PLMD::Value `n`, the forces that PLUMED wishes to apply on the components of the input object are stored in the std::vector called `inputForce`. Furthermore, whenever a PLMD::Value
52 : is created from a PUT action `storedata` is set to true. PUT also has a CONSTANT flag that allows you to transfer variables such as the value of the timestep that is set only once during the
53 : simulation (i.e. during startup).
54 :
55 : Data is transferred from the input void pointers to the PLMD value when the `share` and `wait` methods are called. Vectors, e.g. positions, that are split between the domains
56 : are transferred when the share and wait methods of the [DOMAIN_DECOMPOSITION](DOMAIN_DECOMPOSITION.md) action are called.
57 :
58 : You would only use the PUT command if you were calling PLUMED from python or an MD code. The equivalent commands in a convetional PLUMED input file would look like this.
59 :
60 : ```plumed
61 : # This is how you create a value to hold the energy the MD code passes energy in plumed
62 : eng: PUT UNIT=energy PERIODIC=NO
63 : # This is how you create an vector of the 100 x positions to plumed
64 : # Notice how we use ROLE here in order to tell PLUMED that these are the x coordinates.
65 : # Further note that we need to use the FROM_DOMAINS flag if the MD code we are using uses
66 : # domain decomposition.
67 : xpos: PUT SHAPE=100 UNIT=length PERIODIC=NO ROLE=x FROM_DOMAINS
68 : # This is how you create a scalar to hold the timestep
69 : # The constant flag indicates that the value of the timestep doesn't change during the simulation
70 : tstep: PUT CONSTANT UNIT=time PERIODIC=NO
71 : # This is how you create a value to hold a 10 x 10 matrix in plumed whose elements are unitless
72 : matrix: PUT SHAPE=10,10 UNIT=number PERIODIC=NO
73 : # Lastly, if you want to pass a value that has a periodic domain you can do so as follows
74 : # By adding the MUTABLE flag here we pass the data to PLUMED in a way that ensures that PLUMED can modify
75 : # the value that was passed from the MD code and thus pass back a different value to the underlying MD code.
76 : tor: PUT UNIT=number PERIODIC=-pi,pi MUTABLE
77 : ```
78 :
79 : */
80 : //+ENDPLUMEDOC
81 :
82 : namespace PLMD {
83 :
84 : PLUMED_REGISTER_ACTION(ActionToPutData,"PUT")
85 :
86 13282 : void ActionToPutData::registerKeywords(Keywords& keys) {
87 13282 : ActionForInterface::registerKeywords( keys );
88 13282 : keys.add("compulsory","SHAPE","0","the shape of the value that is being passed to PLUMED");
89 13282 : keys.add("compulsory","UNIT","the unit of the quantity that is being passed to PLUMED through this value. Can be either number, energy, time, length, mass or charge");
90 13282 : keys.add("compulsory","FORCE_UNIT","default","the units to use for the force");
91 13282 : keys.add("compulsory","PERIODIC","if the value being passed to plumed is periodic then you should specify the periodicity of the function. If the value "
92 : "is not periodic you must state this using PERIODIC=NO. Positions are passed with PERIODIC=NO even though special methods are used "
93 : "to deal with pbc");
94 13282 : keys.addFlag("CONSTANT",false,"does this quantity not depend on time");
95 13282 : keys.addFlag("FROM_DOMAINS",false,"is this quantity passed through the domain decomposition object");
96 13282 : keys.addFlag("MUTABLE",false,"can plumed change the value of the pointer that is passed from the MD code");
97 13282 : keys.addFlag("ROUND_TO_BASE_TEN",false,"is this quantity written by the user in base ten, as the timestep is. If the MD code "
98 : "passes floats, the extra digits that appear when converting them to double are removed");
99 13282 : keys.remove("NUMERICAL_DERIVATIVES");
100 26564 : keys.setValueDescription("scalar/vector/matrix/grid","the data that was passed from the MD code");
101 13282 : }
102 :
103 15462 : ActionToPutData::ActionToPutData(const ActionOptions&ao):
104 : Action(ao),
105 : ActionForInterface(ao),
106 15462 : noforce(false),
107 15462 : forces_added_elsewhere(false),
108 15462 : fixed(false),
109 15462 : from_domains(false),
110 15462 : resetable(false),
111 15462 : dataCanBeSet(true),
112 15462 : unit(n),
113 15462 : mydata(DataPassingObject::create(plumed.getRealPrecision())) {
114 15462 : if( getName()!="ENERGY" && getName()!="PBC" ) {
115 13245 : std::vector<std::size_t> shape;
116 26490 : parseVector("SHAPE",shape);
117 13245 : if( shape.size()==1 && shape[0]==0 ) {
118 2327 : shape.resize(0);
119 2327 : addValue( shape );
120 : } else {
121 10918 : addValue( shape );
122 : }
123 :
124 : std::string unitstr, funitstr;
125 13245 : parse("UNIT",unitstr);
126 13245 : parse("FORCE_UNIT",funitstr);
127 13245 : setUnit( unitstr, funitstr );
128 :
129 : // Now sort out period
130 13245 : std::vector<std::string> period;
131 26490 : parseVector("PERIODIC",period);
132 13245 : if( period.size()==1 ) {
133 13245 : if( period[0]!="NO") {
134 0 : error("input to PERIODIC keyword does not make sense");
135 : }
136 13245 : setNotPeriodic();
137 0 : } else if( period.size()==2 ) {
138 0 : setPeriodic( period[0], period[1] );
139 : } else {
140 0 : error("input to PERIODIC keyword does not make sense");
141 : }
142 :
143 13245 : parseFlag("CONSTANT",fixed);
144 13245 : if( fixed ) {
145 6673 : noforce=true;
146 6673 : copyOutput(0)->setConstant();
147 : }
148 13245 : parseFlag("FROM_DOMAINS",from_domains);
149 13245 : parseFlag("MUTABLE",resetable);
150 13245 : bool round_to_base_ten=false;
151 26490 : parseFlag("ROUND_TO_BASE_TEN",round_to_base_ten);
152 : mydata->setRoundToBaseTen( round_to_base_ten );
153 26490 : }
154 15462 : }
155 :
156 15462 : void ActionToPutData::setUnit( const std::string& unitstr, const std::string& funitstr ) {
157 15462 : if( unitstr=="number" ) {
158 52 : unit=n;
159 15410 : } else if( unitstr=="energy" ) {
160 371 : unit=e;
161 15039 : } else if( unitstr=="length" ) {
162 8704 : unit=l;
163 6335 : } else if( unitstr=="mass" ) {
164 2176 : unit=m;
165 4159 : } else if( unitstr=="charge" ) {
166 2176 : unit=q;
167 1983 : } else if( unitstr=="time" ) {
168 1983 : unit=t;
169 : } else {
170 0 : error( unitstr + " is not a valid input unit");
171 : }
172 : // Set the force units
173 15462 : if( funitstr=="default" ) {
174 13286 : funit=d;
175 2176 : } else if( funitstr=="energy" ) {
176 2176 : funit=eng;
177 : } else {
178 0 : error( funitstr + " is not a valid input force unit");
179 : }
180 15462 : }
181 :
182 17580 : std::string ActionToPutData::getUnitName() const {
183 17580 : if( unit==e ) {
184 704 : return "energy";
185 : }
186 : if( unit==l ) {
187 8624 : return "length";
188 : }
189 : if( unit==m ) {
190 2156 : return "mass";
191 : }
192 : if( unit==q ) {
193 2156 : return "charge";
194 : }
195 : if( unit==t ) {
196 3940 : return "time";
197 : }
198 0 : plumed_error();
199 : }
200 :
201 973768 : void ActionToPutData::setStart( const std::string& actname, const unsigned& sss) {
202 973768 : plumed_assert( actname==getLabel() );
203 : mydata->setStart(sss);
204 973768 : }
205 :
206 973768 : void ActionToPutData::setStride( const std::string& actname, const unsigned& sss ) {
207 973768 : plumed_assert( actname==getLabel() );
208 : mydata->setStride(sss);
209 973768 : }
210 :
211 17618 : void ActionToPutData::updateUnits( DataPassingTools* passtools ) {
212 : // Don't need to do anythign if this is just a number
213 17618 : if( unit==n ) {
214 : return ;
215 : }
216 :
217 35160 : double vunits=passtools->getUnitConversion( getUnitName() );
218 : mydata->setUnit(vunits);
219 17580 : if( fixed && wasset ) {
220 4603 : mydata->share_data( 0, getPntrToValue()->getNumberOfValues(), getPntrToValue() );
221 : }
222 17580 : if( funit==eng ) {
223 2156 : mydata->setForceUnit( 1/passtools->getUnitConversion("energy"));
224 15424 : } else if( funit==d ) {
225 15424 : mydata->setForceUnit(1/passtools->getUnitConversion("energy")*vunits);
226 : }
227 : }
228 :
229 3677183 : bool ActionToPutData::setValuePointer( const std::string& actname, const TypesafePtr & val ) {
230 3677183 : if( actname!=getLabel() ) {
231 : return false;
232 : }
233 976087 : wasset=true;
234 976093 : plumed_massert( dataCanBeSet, "set " + getLabel() + " cannot be set at this time");
235 976084 : if( !from_domains ) {
236 221375 : if( !resetable && getPntrToComponent(0)->getRank()==0 ) {
237 6412 : mydata->saveValueAsDouble( val );
238 6412 : if( fixed ) {
239 2323 : mydata->share_data( 0, getPntrToValue()->getNumberOfValues(), getPntrToValue() );
240 : }
241 : } else {
242 214963 : mydata->setValuePointer(val,getPntrToComponent(0)->getShape(), !resetable);
243 : }
244 : } else {
245 754709 : mydata->setValuePointer(val,std::vector<std::size_t>(), !resetable);
246 : }
247 : return true;
248 : }
249 :
250 2011496 : bool ActionToPutData::setForcePointer( const std::string& actname, const TypesafePtr & val ) {
251 2011496 : if( actname!=getLabel() ) {
252 : return false;
253 : }
254 657478 : plumed_massert( dataCanBeSet, "force on " + getLabel() + " cannot be set at this time");
255 657478 : if( !from_domains ) {
256 158182 : mydata->setForcePointer(val,getPntrToComponent(0)->getShape());
257 : } else {
258 499296 : mydata->setForcePointer(val,std::vector<std::size_t>());
259 : }
260 : return true;
261 : }
262 :
263 1422 : void ActionToPutData::getLocalValues( std::vector<double>& vals ) const {
264 1422 : mydata->share_data( vals );
265 1422 : }
266 :
267 219054 : void ActionToPutData::wait() {
268 219054 : dataCanBeSet=false;
269 219054 : if( fixed || !wasset ) {
270 : return;
271 : }
272 : plumed_assert( wasset );
273 216938 : mydata->share_data( 0, getPntrToValue()->getNumberOfValues(), getPntrToValue() );
274 : }
275 :
276 1036729 : void ActionToPutData::apply() {
277 1036729 : if( getPntrToValue()->forcesWereAdded() && !noforce && !forces_added_elsewhere ) {
278 547466 : if( getName()=="ENERGY" || getDependencies().size()==0 ) {
279 129374 : checkForcePointer();
280 129371 : mydata->add_force( getPntrToValue() );
281 : }
282 : }
283 1036726 : }
284 :
285 546415 : void ActionToPutData::checkForcePointer() const {
286 546415 : if( !mydata->hasForcePointer() ) {
287 18 : plumed_merror("PLUMED has to apply forces on " + getLabel() + " but the MD code did not pass the array where these forces should be stored. "
288 : "Pass it with the command setValueForces (setForces for positions, setVirial for the cell) or, if this quantity cannot be biased, make sure no bias acts on it");
289 : }
290 546409 : }
291 :
292 150 : unsigned ActionToPutData::getNumberOfForcesToRescale() const {
293 : // Values that are passed through the domain decomposition only store the atoms
294 : // that are local to this rank, so the number of forces held in the buffer of the
295 : // MD code is the number of local atoms and not the number of values
296 150 : if( !from_domains || getDependencies().size()!=1 ) {
297 0 : return copyOutput(0)->getNumberOfValues();
298 : }
299 150 : plumed_assert(getDependencies()[0]); // needed for following calls, see #1046
300 150 : ActionForInterface* ai = getDependencies()[0]->castToActionForInterface();
301 150 : if( !ai ) {
302 0 : return copyOutput(0)->getNumberOfValues();
303 : }
304 150 : return ai->getNumberOfForcesToRescale();
305 : }
306 :
307 200 : void ActionToPutData::rescaleForces( const double& alpha ) {
308 200 : if( noforce ) {
309 : return;
310 : }
311 150 : wasscaled=true;
312 150 : mydata->rescale_force( getNumberOfForcesToRescale(), alpha, getPntrToValue() );
313 :
314 : }
315 :
316 882 : void ActionToPutData::writeBinary(std::ostream&o) {
317 882 : if(!fixed) {
318 504 : getPntrToValue()->writeBinary(o);
319 : }
320 882 : }
321 :
322 882 : void ActionToPutData::readBinary(std::istream&i) {
323 882 : if(!fixed) {
324 504 : getPntrToValue()->readBinary(i);
325 : }
326 882 : }
327 :
328 : }
|
