ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/applications/staticProps/StaticProps.cpp
Revision: 1999
Committed: Sat May 31 21:00:46 2014 UTC (10 years, 11 months ago) by gezelter
File size: 15741 byte(s)
Log Message:
Updating MultipoleSum to scan the cutoff radii

File Contents

# Content
1 /*
2 * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3 *
4 * The University of Notre Dame grants you ("Licensee") a
5 * non-exclusive, royalty free, license to use, modify and
6 * redistribute this software in source and binary code form, provided
7 * that the following conditions are met:
8 *
9 * 1. Redistributions of source code must retain the above copyright
10 * notice, this list of conditions and the following disclaimer.
11 *
12 * 2. Redistributions in binary form must reproduce the above copyright
13 * notice, this list of conditions and the following disclaimer in the
14 * documentation and/or other materials provided with the
15 * distribution.
16 *
17 * This software is provided "AS IS," without a warranty of any
18 * kind. All express or implied conditions, representations and
19 * warranties, including any implied warranty of merchantability,
20 * fitness for a particular purpose or non-infringement, are hereby
21 * excluded. The University of Notre Dame and its licensors shall not
22 * be liable for any damages suffered by licensee as a result of
23 * using, modifying or distributing the software or its
24 * derivatives. In no event will the University of Notre Dame or its
25 * licensors be liable for any lost revenue, profit or data, or for
26 * direct, indirect, special, consequential, incidental or punitive
27 * damages, however caused and regardless of the theory of liability,
28 * arising out of the use of or inability to use software, even if the
29 * University of Notre Dame has been advised of the possibility of
30 * such damages.
31 *
32 * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your
33 * research, please cite the appropriate papers when you publish your
34 * work. Good starting points are:
35 *
36 * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
37 * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
38 * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).
39 * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010).
40 * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41 */
42
43 #include <iostream>
44 #include <fstream>
45 #include <string>
46
47 #include "brains/SimCreator.hpp"
48 #include "brains/SimInfo.hpp"
49 #include "io/DumpReader.hpp"
50 #include "utils/simError.h"
51
52 #include "applications/staticProps/StaticPropsCmd.h"
53 #include "applications/staticProps/StaticAnalyser.hpp"
54 #include "applications/staticProps/GofR.hpp"
55 #include "applications/staticProps/GofZ.hpp"
56 #include "applications/staticProps/GofRZ.hpp"
57 #include "applications/staticProps/GofRAngle.hpp"
58 #include "applications/staticProps/GofAngle2.hpp"
59 #include "applications/staticProps/GofXyz.hpp"
60 #include "applications/staticProps/TwoDGofR.hpp"
61 #include "applications/staticProps/P2OrderParameter.hpp"
62 #include "applications/staticProps/BondOrderParameter.hpp"
63 #include "applications/staticProps/BOPofR.hpp"
64 #include "applications/staticProps/RippleOP.hpp"
65 #include "applications/staticProps/SCDOrderParameter.hpp"
66 #include "applications/staticProps/DensityPlot.hpp"
67 #include "applications/staticProps/ObjectCount.hpp"
68 #include "applications/staticProps/RhoZ.hpp"
69 #include "applications/staticProps/pAngle.hpp"
70 #include "applications/staticProps/BondAngleDistribution.hpp"
71 #include "applications/staticProps/NanoVolume.hpp"
72 #include "applications/staticProps/NanoLength.hpp"
73 #if defined(HAVE_FFTW_H) || defined(HAVE_DFFTW_H) || defined(HAVE_FFTW3_H)
74 #include "applications/staticProps/Hxy.hpp"
75 #endif
76 #include "applications/staticProps/RhoR.hpp"
77 #include "applications/staticProps/AngleR.hpp"
78 #include "applications/staticProps/TetrahedralityParam.hpp"
79 #include "applications/staticProps/TetrahedralityParamZ.hpp"
80 #include "applications/staticProps/RNEMDStats.hpp"
81 #include "applications/staticProps/NitrileFrequencyMap.hpp"
82 #include "applications/staticProps/MultipoleSum.hpp"
83
84 using namespace OpenMD;
85
86 int main(int argc, char* argv[]){
87
88
89 gengetopt_args_info args_info;
90
91 //parse the command line option
92 if (cmdline_parser (argc, argv, &args_info) != 0) {
93 exit(1) ;
94 }
95
96 //get the dumpfile name
97 std::string dumpFileName = args_info.input_arg;
98 std::string sele1;
99 std::string sele2;
100
101 // check the first selection argument, or set it to the environment
102 // variable, or failing that, set it to "select all"
103
104 if (args_info.sele1_given) {
105 sele1 = args_info.sele1_arg;
106 } else {
107 char* sele1Env= getenv("SELECTION1");
108 if (sele1Env) {
109 sele1 = sele1Env;
110 } else {
111 sele1 = "select all";
112 }
113 }
114
115 // check the second selection argument, or set it to the environment
116 // variable, or failing that, set it to the first selection
117
118 if (args_info.sele2_given) {
119 sele2 = args_info.sele2_arg;
120 } else {
121 char* sele2Env = getenv("SELECTION2");
122 if (sele2Env) {
123 sele2 = sele2Env;
124 } else {
125 //If sele2 is not specified, then the default behavior
126 //should be what is already intended for sele1
127 sele2 = sele1;
128 }
129 }
130
131 bool batchMode;
132 if (args_info.scd_given){
133 if (args_info.sele1_given &&
134 args_info.sele2_given && args_info.sele3_given) {
135 batchMode = false;
136 } else if (args_info.molname_given &&
137 args_info.begin_given && args_info.end_given) {
138 if (args_info.begin_arg < 0 ||
139 args_info.end_arg < 0 || args_info.begin_arg > args_info.end_arg-2) {
140 sprintf( painCave.errMsg,
141 "below conditions are not satisfied:\n"
142 "0 <= begin && 0<= end && begin <= end-2\n");
143 painCave.severity = OPENMD_ERROR;
144 painCave.isFatal = 1;
145 simError();
146 }
147 batchMode = true;
148 } else{
149 sprintf( painCave.errMsg,
150 "either --sele1, --sele2, --sele3 are specified,"
151 " or --molname, --begin, --end are specified\n");
152 painCave.severity = OPENMD_ERROR;
153 painCave.isFatal = 1;
154 simError();
155 }
156 }
157
158 //parse md file and set up the system
159 SimCreator creator;
160 SimInfo* info = creator.createSim(dumpFileName);
161
162 RealType maxLen;
163 RealType zmaxLen;
164 if (args_info.length_given) {
165 maxLen = args_info.length_arg;
166 if (args_info.zlength_given){
167 zmaxLen = args_info.zlength_arg;
168 }
169 } else {
170 Mat3x3d hmat = info->getSnapshotManager()->getCurrentSnapshot()->getHmat();
171 maxLen = std::min(std::min(hmat(0, 0), hmat(1, 1)), hmat(2, 2)) /2.0;
172 zmaxLen = hmat(2,2);
173 }
174
175 StaticAnalyser* analyser;
176 if (args_info.gofr_given){
177 analyser= new GofR(info, dumpFileName, sele1, sele2, maxLen,
178 args_info.nbins_arg);
179 } else if (args_info.gofz_given) {
180 analyser= new GofZ(info, dumpFileName, sele1, sele2, maxLen,
181 args_info.nbins_arg);
182 } else if (args_info.r_z_given) {
183 analyser = new GofRZ(info, dumpFileName, sele1, sele2, maxLen, zmaxLen,
184 args_info.nbins_arg, args_info.nbins_z_arg);
185 } else if (args_info.r_theta_given) {
186 analyser = new GofRTheta(info, dumpFileName, sele1, sele2, maxLen,
187 args_info.nbins_arg, args_info.nanglebins_arg);
188 } else if (args_info.r_omega_given) {
189 analyser = new GofROmega(info, dumpFileName, sele1, sele2, maxLen,
190 args_info.nbins_arg, args_info.nanglebins_arg);
191 } else if (args_info.theta_omega_given) {
192 analyser = new GofAngle2(info, dumpFileName, sele1, sele2,
193 args_info.nanglebins_arg);
194 } else if (args_info.gxyz_given) {
195 if (args_info.refsele_given) {
196 analyser= new GofXyz(info, dumpFileName, sele1, sele2,
197 args_info.refsele_arg, maxLen, args_info.nbins_arg);
198 } else {
199 sprintf( painCave.errMsg,
200 "--refsele must set when --gxyz is used");
201 painCave.severity = OPENMD_ERROR;
202 painCave.isFatal = 1;
203 simError();
204 }
205 } else if (args_info.twodgofr_given){
206 if (args_info.dz_given) {
207 analyser= new TwoDGofR(info, dumpFileName, sele1, sele2, maxLen,
208 args_info.dz_arg, args_info.nbins_arg);
209 } else {
210 sprintf( painCave.errMsg,
211 "A slab width (dz) must be specified when calculating TwoDGofR");
212 painCave.severity = OPENMD_ERROR;
213 painCave.isFatal = 1;
214 simError();
215 }
216 } else if (args_info.p2_given) {
217 if (args_info.sele1_given) {
218 if (args_info.sele2_given)
219 analyser = new P2OrderParameter(info, dumpFileName, sele1, sele2);
220 else
221 if (args_info.seleoffset_given)
222 analyser = new P2OrderParameter(info, dumpFileName, sele1,
223 args_info.seleoffset_arg);
224 else
225 analyser = new P2OrderParameter(info, dumpFileName, sele1);
226 } else {
227 sprintf( painCave.errMsg,
228 "At least one selection script (--sele1) must be specified when calculating P2 order parameters");
229 painCave.severity = OPENMD_ERROR;
230 painCave.isFatal = 1;
231 simError();
232 }
233 } else if (args_info.rp2_given){
234 analyser = new RippleOP(info, dumpFileName, sele1, sele2);
235 } else if (args_info.bo_given){
236 if (args_info.rcut_given) {
237 analyser = new BondOrderParameter(info, dumpFileName, sele1,
238 args_info.rcut_arg,
239 args_info.nbins_arg);
240 } else {
241 sprintf( painCave.errMsg,
242 "A cutoff radius (rcut) must be specified when calculating Bond Order Parameters");
243 painCave.severity = OPENMD_ERROR;
244 painCave.isFatal = 1;
245 simError();
246 }
247 } else if (args_info.multipole_given){
248 analyser = new MultipoleSum(info, dumpFileName, sele1,
249 maxLen, args_info.nbins_arg);
250 } else if (args_info.tet_param_given) {
251 if (args_info.rcut_given) {
252 analyser = new TetrahedralityParam(info, dumpFileName, sele1,
253 args_info.rcut_arg,
254 args_info.nbins_arg);
255 } else {
256 sprintf( painCave.errMsg,
257 "A cutoff radius (rcut) must be specified when calculating Tetrahedrality Parameters");
258 painCave.severity = OPENMD_ERROR;
259 painCave.isFatal = 1;
260 simError();
261 }
262 } else if (args_info.tet_param_z_given) {
263 if (args_info.rcut_given) {
264 analyser = new TetrahedralityParamZ(info, dumpFileName, sele1, sele2,
265 args_info.rcut_arg,
266 args_info.nbins_arg);
267 } else {
268 sprintf( painCave.errMsg,
269 "A cutoff radius (rcut) must be specified when calculating Tetrahedrality Parameters");
270 painCave.severity = OPENMD_ERROR;
271 painCave.isFatal = 1;
272 simError();
273 }
274 } else if (args_info.ior_given){
275 if (args_info.rcut_given) {
276 analyser = new IcosahedralOfR(info, dumpFileName, sele1,
277 args_info.rcut_arg,
278 args_info.nbins_arg, maxLen);
279 } else {
280 sprintf( painCave.errMsg,
281 "A cutoff radius (rcut) must be specified when calculating Bond Order Parameters");
282 painCave.severity = OPENMD_ERROR;
283 painCave.isFatal = 1;
284 simError();
285 }
286 } else if (args_info.for_given){
287 if (args_info.rcut_given) {
288 analyser = new FCCOfR(info, dumpFileName, sele1, args_info.rcut_arg,
289 args_info.nbins_arg, maxLen);
290 } else {
291 sprintf( painCave.errMsg,
292 "A cutoff radius (rcut) must be specified when calculating Bond Order Parameters");
293 painCave.severity = OPENMD_ERROR;
294 painCave.isFatal = 1;
295 simError();
296 }
297 } else if (args_info.bad_given){
298 if (args_info.rcut_given) {
299 analyser = new BondAngleDistribution(info, dumpFileName, sele1,
300 args_info.rcut_arg,
301 args_info.nbins_arg);
302 } else {
303 sprintf( painCave.errMsg,
304 "A cutoff radius (rcut) must be specified when calculating Bond Angle Distributions");
305 painCave.severity = OPENMD_ERROR;
306 painCave.isFatal = 1;
307 simError();
308 }
309 } else if (args_info.scd_given) {
310 if (batchMode) {
311 analyser = new SCDOrderParameter(info, dumpFileName,
312 args_info.molname_arg,
313 args_info.begin_arg, args_info.end_arg);
314 } else{
315 std::string sele3 = args_info.sele3_arg;
316 analyser = new SCDOrderParameter(info, dumpFileName,
317 sele1, sele2, sele3);
318 }
319 }else if (args_info.density_given) {
320 analyser= new DensityPlot(info, dumpFileName, sele1, sele2, maxLen,
321 args_info.nbins_arg);
322 } else if (args_info.count_given) {
323 analyser = new ObjectCount(info, dumpFileName, sele1 );
324 } else if (args_info.slab_density_given) {
325 analyser = new RhoZ(info, dumpFileName, sele1, args_info.nbins_arg);
326 } else if (args_info.rnemdz_given) {
327 analyser = new RNEMDZ(info, dumpFileName, sele1, args_info.nbins_arg);
328 } else if (args_info.rnemdr_given) {
329 analyser = new RNEMDR(info, dumpFileName, sele1, args_info.nbins_arg);
330 } else if (args_info.rnemdrt_given) {
331 analyser = new RNEMDRTheta(info, dumpFileName, sele1,
332 args_info.nbins_arg, args_info.nanglebins_arg);
333 } else if (args_info.nitrile_given) {
334 analyser = new NitrileFrequencyMap(info, dumpFileName, sele1,
335 args_info.nbins_arg);
336 } else if (args_info.p_angle_given) {
337 if (args_info.sele1_given) {
338 if (args_info.sele2_given)
339 analyser = new pAngle(info, dumpFileName, sele1, sele2,
340 args_info.nbins_arg);
341 else
342 if (args_info.seleoffset_given) {
343 if (args_info.seleoffset2_given) {
344 analyser = new pAngle(info, dumpFileName, sele1,
345 args_info.seleoffset_arg,
346 args_info.seleoffset2_arg,
347 args_info.nbins_arg);
348 } else {
349 analyser = new pAngle(info, dumpFileName, sele1,
350 args_info.seleoffset_arg,
351 args_info.nbins_arg);
352 }
353 } else
354 analyser = new pAngle(info, dumpFileName, sele1,
355 args_info.nbins_arg);
356 } else {
357 sprintf( painCave.errMsg,
358 "At least one selection script (--sele1) must be specified when "
359 "calculating P(angle) distributions");
360 painCave.severity = OPENMD_ERROR;
361 painCave.isFatal = 1;
362 simError();
363 }
364 #if defined(HAVE_FFTW_H) || defined(HAVE_DFFTW_H) || defined(HAVE_FFTW3_H)
365 }else if (args_info.hxy_given) {
366 analyser = new Hxy(info, dumpFileName, sele1, args_info.nbins_x_arg,
367 args_info.nbins_y_arg, args_info.nbins_arg);
368 #endif
369 }else if (args_info.rho_r_given) {
370 if (args_info.radius_given){
371 analyser = new RhoR(info, dumpFileName, sele1, maxLen,args_info.nbins_arg,args_info.radius_arg);
372 }else{
373 sprintf( painCave.errMsg,
374 "A particle radius (radius) must be specified when calculating Rho(r)");
375 painCave.severity = OPENMD_ERROR;
376 painCave.isFatal = 1;
377 simError();
378 }
379 } else if (args_info.hullvol_given) {
380 analyser = new NanoVolume(info, dumpFileName, sele1);
381 } else if (args_info.rodlength_given) {
382 analyser = new NanoLength(info, dumpFileName, sele1);
383 } else if (args_info.angle_r_given) {
384 analyser = new AngleR(info, dumpFileName, sele1, maxLen,args_info.nbins_arg);
385 }
386
387 if (args_info.output_given) {
388 analyser->setOutputName(args_info.output_arg);
389 }
390 if (args_info.step_given) {
391 analyser->setStep(args_info.step_arg);
392 }
393
394 analyser->process();
395
396 delete analyser;
397 delete info;
398
399 return 0;
400 }

Properties

Name Value
svn:executable *
svn:keywords Author Id Revision Date