OpenMD 3.0
Molecular Dynamics in the Open
Loading...
Searching...
No Matches
DynamicProps.cpp
1/*
2 * Copyright (c) 2004-present, The University of Notre Dame. All rights
3 * reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * 1. Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 *
11 * 2. Redistributions in binary form must reproduce the above copyright notice,
12 * this list of conditions and the following disclaimer in the documentation
13 * and/or other materials provided with the distribution.
14 *
15 * 3. Neither the name of the copyright holder nor the names of its
16 * contributors may be used to endorse or promote products derived from
17 * this software without specific prior written permission.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
23 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
24 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
25 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
26 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
27 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
28 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
29 * POSSIBILITY OF SUCH DAMAGE.
30 *
31 * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your
32 * research, please cite the appropriate papers when you publish your
33 * work. Good starting points are:
34 *
35 * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
36 * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
37 * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).
38 * [4] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
39 * [5] Kuang & Gezelter, Mol. Phys., 110, 691-701 (2012).
40 * [6] Lamichhane, Gezelter & Newman, J. Chem. Phys. 141, 134109 (2014).
41 * [7] Lamichhane, Newman & Gezelter, J. Chem. Phys. 141, 134110 (2014).
42 * [8] Bhattarai, Newman & Gezelter, Phys. Rev. B 99, 094106 (2019).
43 */
44
45#include <fstream>
46#include <iostream>
47#include <memory>
48#include <string>
49
50#include "DynamicPropsCmd.hpp"
51#include "applications/dynamicProps/AngularVelVelOutProdCorrFunc.hpp"
52#include "applications/dynamicProps/AngularVelocityAutoOutProductCorrFunc.hpp"
53#include "applications/dynamicProps/BondCorrFunc.hpp"
54#include "applications/dynamicProps/ChargeKineticCorrFunc.hpp"
55#include "applications/dynamicProps/ChargeOrientationCorrFunc.hpp"
56#include "applications/dynamicProps/CollectiveDipoleDisplacement.hpp"
57#include "applications/dynamicProps/CurrentDensityAutoCorrFunc.hpp"
58#include "applications/dynamicProps/DipoleCorrFunc.hpp"
59#include "applications/dynamicProps/DirectionalRCorrFunc.hpp"
60#include "applications/dynamicProps/Displacement.hpp"
61#include "applications/dynamicProps/ForTorCorrFunc.hpp"
62#include "applications/dynamicProps/ForceAutoCorrFunc.hpp"
64#include "applications/dynamicProps/HBondJump.hpp"
65#include "applications/dynamicProps/HBondPersistence.hpp"
66#include "applications/dynamicProps/LegendreCorrFunc.hpp"
67#include "applications/dynamicProps/LegendreCorrFuncZ.hpp"
68#include "applications/dynamicProps/MomAngMomCorrFunc.hpp"
69#include "applications/dynamicProps/OnsagerCorrFunc.hpp"
70#include "applications/dynamicProps/RCorrFunc.hpp"
71#include "applications/dynamicProps/RotAngleDisplacement.hpp"
72#include "applications/dynamicProps/SelectionCorrFunc.hpp"
73#include "applications/dynamicProps/StressCorrFunc.hpp"
74#include "applications/dynamicProps/SystemDipoleCorrFunc.hpp"
75#include "applications/dynamicProps/ThetaCorrFunc.hpp"
76#include "applications/dynamicProps/TorForCorrFunc.hpp"
77#include "applications/dynamicProps/TorqueAutoCorrFunc.hpp"
78#include "applications/dynamicProps/VCorrFunc.hpp"
79#include "applications/dynamicProps/VelAngularVelOutProdCorrFunc.hpp"
80#include "applications/dynamicProps/VelocityAutoOutProductCorrFunc.hpp"
81#include "applications/dynamicProps/WCorrFunc.hpp"
82#include "applications/dynamicProps/cOHz.hpp"
83#include "brains/SimCreator.hpp"
84#include "brains/SimInfo.hpp"
85#include "utils/Revision.hpp"
86#include "utils/StringUtils.hpp"
87#include "utils/simError.h"
88
89using namespace OpenMD;
90
91int main(int argc, char* argv[]) {
92 struct gengetopt_args_info args_info;
93
94 // parse the command line option
95
96 if (cmdline_parser(argc, argv, &args_info) != 0) {
98 exit(1);
99 }
100
101 // get the dumpfile name and meta-data file name
102 std::string dumpFileName = args_info.input_arg;
103
104 std::string sele1;
105 std::string sele2;
106
107 // check the first selection argument, or set it to the environment
108 // variable, or failing that, set it to "select all"
109
110 if (args_info.sele1_given) {
111 sele1 = args_info.sele1_arg;
112 } else {
113 char* sele1Env = getenv("SELECTION1");
114 if (sele1Env) {
115 sele1 = sele1Env;
116 } else {
117 sele1 = "select all";
118 }
119 }
120
121 // check the second selection argument, or set it to the environment
122 // variable, or failing that, set it to the first selection
123
124 if (args_info.sele2_given) {
125 sele2 = args_info.sele2_arg;
126 } else {
127 char* sele2Env = getenv("SELECTION2");
128 if (sele2Env) {
129 sele2 = sele2Env;
130 } else {
131 // If sele2 is not specified, then the default behavior
132 // should be what is already intended for sele1
133 sele2 = sele1;
134 }
135 }
136
137 // convert privilegedAxis to corresponding integer
138 // x axis -> 0
139 // y axis -> 1
140 // z axis -> 2 (default)
141
142 int privilegedAxis;
143 switch (args_info.privilegedAxis_arg) {
144 case privilegedAxis_arg_x:
145 privilegedAxis = 0;
146 break;
147 case privilegedAxis_arg_y:
148 privilegedAxis = 1;
149 break;
150 case privilegedAxis_arg_z:
151 default:
152 privilegedAxis = 2;
153 break;
154 }
155
156 // use the memory string to figure out how much memory we can use:
157 // char *end;
158 // long long int memSize = memparse(args_info.memory_arg, &end);
159
160 // We don't really need to print this out anymore:
161 // snprintf( painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
162 // "Amount of memory being used: %llu bytes\n", memSize);
163 // painCave.severity = OPENMD_INFO;
164 // painCave.isFatal = 0;
165 // simError();
166
167 // parse md file and set up the system
168 SimCreator creator;
169 SimInfo* info = creator.createSim(dumpFileName, false);
170
171 RealType maxLen;
172 if (args_info.length_given) {
173 maxLen = args_info.length_arg;
174 } else {
175 maxLen = 100.0;
176 }
177
178 std::unique_ptr<DynamicProperty> corrFunc {nullptr};
179
180 if (args_info.sdcorr_given) {
181 corrFunc = std::make_unique<SystemDipoleCorrFunc>(info, dumpFileName, sele1,
182 sele2);
183 } else if (args_info.selecorr_given) {
184 corrFunc =
185 std::make_unique<SelectionCorrFunc>(info, dumpFileName, sele1, sele2);
186 } else if (args_info.dcorr_given) {
187 corrFunc =
188 std::make_unique<DipoleCorrFunc>(info, dumpFileName, sele1, sele2);
189 } else if (args_info.rcorr_given) {
190 corrFunc = std::make_unique<RCorrFunc>(info, dumpFileName, sele1, sele2);
191 } else if (args_info.r_rcorr_given) {
192 corrFunc = std::make_unique<RCorrFuncR>(info, dumpFileName, sele1, sele2);
193 } else if (args_info.thetacorr_given) {
194 corrFunc =
195 std::make_unique<ThetaCorrFunc>(info, dumpFileName, sele1, sele2);
196 } else if (args_info.drcorr_given) {
197 corrFunc = std::make_unique<DirectionalRCorrFunc>(info, dumpFileName, sele1,
198 sele2);
199 } else if (args_info.rcorrZ_given) {
200 corrFunc = std::make_unique<RCorrFuncZ>(
201 info, dumpFileName, sele1, sele2, args_info.nzbins_arg, privilegedAxis);
202 } else if (args_info.vcorr_given) {
203 corrFunc = std::make_unique<VCorrFunc>(info, dumpFileName, sele1, sele2);
204 } else if (args_info.vcorrZ_given) {
205 corrFunc = std::make_unique<VCorrFuncZ>(info, dumpFileName, sele1, sele2);
206 } else if (args_info.vcorrR_given) {
207 corrFunc = std::make_unique<VCorrFuncR>(info, dumpFileName, sele1, sele2);
208 } else if (args_info.wcorr_given) {
209 corrFunc = std::make_unique<WCorrFunc>(info, dumpFileName, sele1, sele2);
210 } else if (args_info.pjcorr_given) {
211 corrFunc =
212 std::make_unique<MomAngMomCorrFunc>(info, dumpFileName, sele1, sele2);
213 } else if (args_info.ftcorr_given) {
214 corrFunc =
215 std::make_unique<ForTorCorrFunc>(info, dumpFileName, sele1, sele2);
216 } else if (args_info.ckcorr_given) {
217 corrFunc = std::make_unique<ChargeKineticCorrFunc>(
218 info, dumpFileName, sele1, sele2, args_info.rcut_arg);
219 } else if (args_info.cscorr_given) {
220 if (args_info.dipoleX_given && args_info.dipoleY_given &&
221 args_info.dipoleZ_given) {
222 corrFunc = std::make_unique<ChargeOrientationCorrFunc>(
223 info, dumpFileName, sele1, sele2, args_info.dipoleX_arg,
224 args_info.dipoleY_arg, args_info.dipoleZ_arg, args_info.rcut_arg);
225 }
226 } else if (args_info.facorr_given) {
227 corrFunc =
228 std::make_unique<ForceAutoCorrFunc>(info, dumpFileName, sele1, sele2);
229 } else if (args_info.tfcorr_given) {
230 corrFunc =
231 std::make_unique<TorForCorrFunc>(info, dumpFileName, sele1, sele2);
232 } else if (args_info.tacorr_given) {
233 corrFunc =
234 std::make_unique<TorqueAutoCorrFunc>(info, dumpFileName, sele1, sele2);
235 } else if (args_info.bondcorr_given) {
236 corrFunc = std::make_unique<BondCorrFunc>(info, dumpFileName, sele1, sele2);
237 } else if (args_info.stresscorr_given) {
238 corrFunc =
239 std::make_unique<StressCorrFunc>(info, dumpFileName, sele1, sele2);
240 } else if (args_info.freqfluccorr_given) {
241 corrFunc =
242 std::make_unique<FreqFlucCorrFunc>(info, dumpFileName, sele1, sele2);
243 } else if (args_info.lcorr_given) {
244 int order(0);
245 if (args_info.order_given)
246 order = args_info.order_arg;
247 else {
248 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
249 "--order must be set if --lcorr is set\n");
250 painCave.severity = OPENMD_ERROR;
251 painCave.isFatal = 1;
252 simError();
253 }
254
255 corrFunc = std::make_unique<LegendreCorrFunc>(info, dumpFileName, sele1,
256 sele2, order);
257 } else if (args_info.lcorrZ_given) {
258 int order(0);
259 if (args_info.order_given)
260 order = args_info.order_arg;
261 else {
262 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
263 "--order must be set if --lcorrZ is set\n");
264 painCave.severity = OPENMD_ERROR;
265 painCave.isFatal = 1;
266 simError();
267 }
268
269 corrFunc = std::make_unique<LegendreCorrFuncZ>(
270 info, dumpFileName, sele1, sele2, order, args_info.nzbins_arg,
271 privilegedAxis);
272
273 } else if (args_info.cohZ_given) {
274 int order(0);
275 if (args_info.order_given)
276 order = args_info.order_arg;
277 else {
278 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
279 "--order must be set if --cohZ is set\n");
280 painCave.severity = OPENMD_ERROR;
281 painCave.isFatal = 1;
282 simError();
283 }
284
285 corrFunc = std::make_unique<COHZ>(info, dumpFileName, sele1, sele2, order,
286 args_info.nzbins_arg, privilegedAxis);
287
288 } else if (args_info.jumptime_given) {
289 corrFunc = std::make_unique<HBondJump>(
290 info, dumpFileName, sele1, sele2, args_info.OOcut_arg,
291 args_info.thetacut_arg, args_info.OHcut_arg);
292 } else if (args_info.jumptimeZ_given) {
293 corrFunc = std::make_unique<HBondJumpZ>(
294 info, dumpFileName, sele1, sele2, args_info.OOcut_arg,
295 args_info.thetacut_arg, args_info.OHcut_arg, args_info.nzbins_arg,
296 privilegedAxis);
297 } else if (args_info.jumptimeR_given) {
298 if (args_info.sele3_given) {
299 corrFunc = std::make_unique<HBondJumpR>(
300 info, dumpFileName, sele1, sele2, args_info.sele3_arg,
301 args_info.OOcut_arg, args_info.thetacut_arg, args_info.OHcut_arg,
302 maxLen, args_info.nbins_arg);
303 } else {
304 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
305 "--sele3 must be set if --jumptimeR is set\n");
306 painCave.severity = OPENMD_ERROR;
307 painCave.isFatal = 1;
308 simError();
309 }
310 } else if (args_info.persistence_given) {
311 corrFunc = std::make_unique<HBondPersistence>(
312 info, dumpFileName, sele1, sele2, args_info.OOcut_arg,
313 args_info.thetacut_arg, args_info.OHcut_arg);
314 } else if (args_info.disp_given) {
315 corrFunc = std::make_unique<Displacement>(info, dumpFileName, sele1, sele2);
316 } else if (args_info.dispZ_given) {
317 corrFunc = std::make_unique<DisplacementZ>(
318 info, dumpFileName, sele1, sele2, args_info.nzbins_arg, privilegedAxis);
319 } else if (args_info.current_given) {
320 corrFunc = std::make_unique<CurrentDensityAutoCorrFunc>(info, dumpFileName,
321 sele1, sele2);
322 } else if (args_info.onsager_given) {
323 if (args_info.sele1_given) {
324 corrFunc =
325 std::make_unique<OnsagerCorrFunc>(info, dumpFileName, sele1, sele2);
326 } else {
327 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
328 "--sele1 must be set for Center of Mass Rcorr\n");
329 painCave.severity = OPENMD_ERROR;
330 painCave.isFatal = 1;
331 simError();
332 }
333 } else if (args_info.ddisp_given) {
334 corrFunc = std::make_unique<CollectiveDipoleDisplacement>(
335 info, dumpFileName, sele1, sele2);
336 } else if (args_info.vaOutProdcorr_given) {
337 corrFunc = std::make_unique<VelocityAutoOutProductCorrFunc>(
338 info, dumpFileName, sele1, sele2);
339 } else if (args_info.waOutProdcorr_given) {
340 corrFunc = std::make_unique<AngularVelocityAutoOutProductCorrFunc>(
341 info, dumpFileName, sele1, sele2);
342 } else if (args_info.vwOutProdcorr_given) {
343 corrFunc = std::make_unique<VelAngularVelOutProdCorrFunc>(
344 info, dumpFileName, sele1, sele2);
345 } else if (args_info.wvOutProdcorr_given) {
346 corrFunc = std::make_unique<AngularVelVelOutProdCorrFunc>(
347 info, dumpFileName, sele1, sele2);
348 } else if (args_info.rotAngleDisp_given) {
349 corrFunc = std::make_unique<RotAngleDisplacement>(info, dumpFileName, sele1,
350 sele2);
351 }
352
353 if (args_info.output_given) { corrFunc->setOutputName(args_info.output_arg); }
354
355 corrFunc->doCorrelate();
356
357 delete info;
358 return 0;
359}
The header file for the command line option parser generated by GNU Gengetopt version 2....
Frequency Fluctuation Correlation Function.
The only responsibility of SimCreator is to parse the meta-data file and create a SimInfo instance ba...
One of the heavy-weight classes of OpenMD, SimInfo maintains objects and variables relating to the cu...
Definition SimInfo.hpp:93
This basic Periodic Table class was originally taken from the data.cpp file in OpenBabel.
Where the command line options are stored.
int cmdline_parser(int argc, char **argv, struct gengetopt_args_info *args_info)
The command line parser.
void cmdline_parser_print_help(void)
Print the help.