ViewVC Help

View File | Revision Log | Show Annotations | View Changeset | Root Listing

root/OpenMD/trunk/src/applications/nanoparticleBuilder/nanoparticleBuilder.cpp

Comparing trunk/src/applications/nanoparticleBuilder/nanoparticleBuilder.cpp (file contents):
Revision 1075 by gezelter, Tue Oct 17 17:51:52 2006 UTC vs.
Revision 1977 by gezelter, Wed Mar 12 21:35:23 2014 UTC

#	Line 6 | Line 6
6		* redistribute this software in source and binary code form, provided
7		* that the following conditions are met:
8		*
9	<	* 1. Acknowledgement of the program authors must be made in any
10	<	*    publication of scientific results based in part on use of the
11	<	*    program.  An acceptable form of acknowledgement is citation of
12	<	*    the article in which the program was described (Matthew
13	<	*    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14	<	*    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15	<	*    Parallel Simulation Engine for Molecular Dynamics,"
16	<	*    J. Comput. Chem. 26, pp. 252-271 (2005))
17	<	*
18	<	* 2. Redistributions of source code must retain the above copyright
9	>	* 1. Redistributions of source code must retain the above copyright
10		*    notice, this list of conditions and the following disclaimer.
11		*
12	<	* 3. Redistributions in binary form must reproduce the above copyright
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.
#	Line 37 | Line 28
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 <cstdlib>
#	Line 64 | Line 65 | using namespace std;
65		#include "utils/StringUtils.hpp"
66
67		using namespace std;
68	<	using namespace oopse;
68	>	using namespace OpenMD;
69		void createMdFile(const std::string&oldMdFileName,
70		const std::string&newMdFileName,
71		std::vector<int> numMol);
72
73		int main(int argc, char *argv []) {
74
74	–	//register force fields
75	–	registerForceFields();
75		registerLattice();
76
77		gengetopt_args_info args_info;
78		std::string latticeType;
79		std::string inputFileName;
80		std::string outputFileName;
82	–
81		MoLocator* locator;
82		int nComponents;
83		double latticeConstant;
84	<	std::vector<double> lc;
87	<
88	<	double particleRadius;
89	<
84	>	RealType particleRadius;
85		Mat3x3d hmat;
91	–	std::vector<Vector3d> latticePos;
92	–	std::vector<Vector3d> latticeOrt;
93	–
86		DumpWriter *writer;
87
88		// Parse Command Line Arguments
#	Line 104 | Line 96 | int main(int argc, char *argv []) {
96		if (args_info.inputs_num)
97		inputFileName = args_info.inputs[0];
98		else {
99	<	sprintf(painCave.errMsg, "No input .md file name was specified"
99	>	sprintf(painCave.errMsg, "No input .md file name was specified "
100		"on the command line");
101		painCave.isFatal = 1;
102		cmdline_parser_print_help();
#	Line 115 | Line 107 | int main(int argc, char *argv []) {
107		SimCreator oldCreator;
108		SimInfo* oldInfo = oldCreator.createSim(inputFileName, false);
109
110	<	latticeConstant = args_info.latticeCnst_arg;
110	>	latticeConstant = args_info.latticeConstant_arg;
111		particleRadius = args_info.radius_arg;
112		Globals* simParams = oldInfo->getSimParams();
113
#	Line 127 | Line 119 | int main(int argc, char *argv []) {
119		vector<Vector3d> sites = nanoParticle.getSites();
120		vector<Vector3d> orientations = nanoParticle.getOrientations();
121
130	–	std::cout <<"nSites: " << sites.size() << std::endl;
122
123	+	std::vector<int> vacancyTargets;
124	+	vector<bool> isVacancy;
125	+
126	+	Vector3d myLoc;
127	+	RealType myR;
128	+
129	+	for (unsigned int i = 0; i < sites.size(); i++)
130	+	isVacancy.push_back(false);
131	+
132	+	if (args_info.vacancyPercent_given) {
133	+	if (args_info.vacancyPercent_arg < 0.0 || args_info.vacancyPercent_arg > 100.0) {
134	+	sprintf(painCave.errMsg, "vacancyPercent was set to a non-sensical value.");
135	+	painCave.isFatal = 1;
136	+	simError();
137	+	} else {
138	+	RealType vF = args_info.vacancyPercent_arg / 100.0;
139	+	RealType vIR;
140	+	RealType vOR;
141	+	if (args_info.vacancyInnerRadius_given) {
142	+	vIR = args_info.vacancyInnerRadius_arg;
143	+	} else {
144	+	vIR = 0.0;
145	+	}
146	+	if (args_info.vacancyOuterRadius_given) {
147	+	vOR = args_info.vacancyOuterRadius_arg;
148	+	} else {
149	+	vOR = particleRadius;
150	+	}
151	+	if (vIR >= 0.0 && vOR <= particleRadius && vOR >= vIR) {
152	+
153	+	for (unsigned int i = 0; i < sites.size(); i++) {
154	+	myLoc = sites[i];
155	+	myR = myLoc.length();
156	+	if (myR >= vIR && myR <= vOR) {
157	+	vacancyTargets.push_back(i);
158	+	}
159	+	}
160	+	std::random_shuffle(vacancyTargets.begin(), vacancyTargets.end());
161	+
162	+	int nTargets = vacancyTargets.size();
163	+	vacancyTargets.resize((int)(vF * nTargets));
164	+
165	+
166	+	sprintf(painCave.errMsg, "Removing %d atoms from randomly-selected\n"
167	+	"\tsites between %lf and %lf.", (int) vacancyTargets.size(),
168	+	vIR, vOR);
169	+	painCave.isFatal = 0;
170	+	painCave.severity = OPENMD_INFO;
171	+	simError();
172	+
173	+	isVacancy.clear();
174	+	for (unsigned int i = 0; i < sites.size(); i++) {
175	+	bool vac = false;
176	+	for (unsigned int j = 0; j < vacancyTargets.size(); j++) {
177	+	if (i == vacancyTargets[j]) vac = true;
178	+	}
179	+	isVacancy.push_back(vac);
180	+	}
181	+
182	+	} else {
183	+	sprintf(painCave.errMsg, "Something is strange about the vacancy\n"
184	+	"\tinner or outer radii.  Check their values.");
185	+	painCave.isFatal = 1;
186	+	simError();
187	+	}
188	+	}
189	+	}
190	+
191		/* Get number of lattice sites */
192	<	int nSites = sites.size();
192	>	int nSites = sites.size() - vacancyTargets.size();
193
194		std::vector<Component*> components = simParams->getComponents();
195		std::vector<RealType> molFractions;
196		std::vector<RealType> shellRadii;
138	–	std::vector<RealType> molecularMasses;
197		std::vector<int> nMol;
198		std::map<int, int> componentFromSite;
199		nComponents = components.size();
200
201	<	if (args_info.molFraction_given && args_info.ShellRadius_given) {
202	<	sprintf(painCave.errMsg, "Specify either molFraction or ShellRadius "
201	>	if (args_info.molFraction_given && args_info.shellRadius_given) {
202	>	sprintf(painCave.errMsg, "Specify either molFraction or shellRadius "
203		"arguments, but not both!");
204		painCave.isFatal = 1;
205		simError();
#	Line 168 | Line 226 | int main(int argc, char *argv []) {
226		painCave.isFatal = 1;
227		simError();
228		}
229	<	} else if ((int)args_info.ShellRadius_given) {
230	<	if ((int)args_info.ShellRadius_given == nComponents) {
229	>	} else if ((int)args_info.shellRadius_given) {
230	>	if ((int)args_info.shellRadius_given == nComponents) {
231		for (int i = 0; i < nComponents; i++) {
232	<	shellRadii.push_back(args_info.ShellRadius_arg[i]);
232	>	shellRadii.push_back(args_info.shellRadius_arg[i]);
233		}
234	<	} else if ((int)args_info.ShellRadius_given == nComponents-1) {
234	>	} else if ((int)args_info.shellRadius_given == nComponents-1) {
235		for (int i = 0; i < nComponents-1; i++) {
236	<	shellRadii.push_back(args_info.ShellRadius_arg[i]);
236	>	shellRadii.push_back(args_info.shellRadius_arg[i]);
237		}
238		shellRadii.push_back(particleRadius);
239		} else {
240	<	sprintf(painCave.errMsg, "nanoparticleBuilder can't figure out the shell radii "
241	<	"for all of the components in the <MetaData> block.");
240	>	sprintf(painCave.errMsg, "nanoparticleBuilder can't figure out the\n"
241	>	"\tshell radii for all of the components in the <MetaData> block.");
242		painCave.isFatal = 1;
243		simError();
244		}
245		} else {
246	<	sprintf(painCave.errMsg, "You have a multi-component <MetaData> block, but have not "
247	<	"specified either molFraction or ShellRadius arguments.");
246	>	sprintf(painCave.errMsg, "You have a multi-component <MetaData> block,\n"
247	>	"\tbut have not specified either molFraction or shellRadius arguments.");
248		painCave.isFatal = 1;
249		simError();
250		}
#	Line 240 | Line 298 | int main(int argc, char *argv []) {
298		}
299		} else {
300
301	<	for (int i = 0; i < shellRadii.size(); i++) {
301	>	for (unsigned int i = 0; i < shellRadii.size(); i++) {
302		if (shellRadii.at(i) > particleRadius + 1e-6 ) {
303		sprintf(painCave.errMsg, "One of the shellRadius values exceeds the particle Radius.");
304		painCave.isFatal = 1;
#	Line 253 | Line 311 | int main(int argc, char *argv []) {
311		}
312		}
313		}
314	<
315	<	vector<int> ids;
258	<	for (int i = 0; i < sites.size(); i++) ids.push_back(i);
259	<	/* Random particle is the default case*/
314	>
315	>	vector<int> ids;
316		if ((int)args_info.molFraction_given){
317		sprintf(painCave.errMsg, "Creating a randomized spherical nanoparticle.");
318		painCave.isFatal = 0;
319	+	painCave.severity = OPENMD_INFO;
320		simError();
321	+	/* Random particle is the default case*/
322	+
323	+	for (unsigned int i = 0; i < sites.size(); i++)
324	+	if (!isVacancy[i]) ids.push_back(i);
325	+
326		std::random_shuffle(ids.begin(), ids.end());
327	+
328		} else{
329		sprintf(painCave.errMsg, "Creating a core-shell spherical nanoparticle.");
330		painCave.isFatal = 0;
331	+	painCave.severity = OPENMD_INFO;
332		simError();
333
270	–	Vector3d myLoc;
271	–	RealType myR;
334		RealType smallestSoFar;
335		int myComponent = -1;
336		nMol.clear();
337		nMol.resize(nComponents);
338
339	<	for (int i = 0; i < sites.size(); i++) {
339	>	for (unsigned int i = 0; i < sites.size(); i++) {
340		myLoc = sites[i];
341		myR = myLoc.length();
342	<	smallestSoFar = particleRadius;
343	<
344	<	for (int j = 0; j < nComponents; j++) {
345	<	if (myR <= shellRadii[j]) {
346	<	if (shellRadii[j] <= smallestSoFar) {
347	<	smallestSoFar = shellRadii[j];
348	<	myComponent = j;
342	>	smallestSoFar = particleRadius;
343	>	if (!isVacancy[i]) {
344	>	for (int j = 0; j < nComponents; j++) {
345	>	if (myR <= shellRadii[j]) {
346	>	if (shellRadii[j] <= smallestSoFar) {
347	>	smallestSoFar = shellRadii[j];
348	>	myComponent = j;
349	>	}
350		}
351		}
352	+	componentFromSite[i] = myComponent;
353	+	nMol[myComponent]++;
354		}
355	<	componentFromSite[i] = myComponent;
356	<	nMol[myComponent]++;
292	<	}
293	<	}
355	>	}
356	>	}
357
358		outputFileName = args_info.output_arg;
359	<
297	<
359	>
360		//creat new .md file on fly which corrects the number of molecule
361		createMdFile(inputFileName, outputFileName, nMol);
362
363	<	if (oldInfo != NULL)
302	<	delete oldInfo;
363	>	delete oldInfo;
364
365		SimCreator newCreator;
366		SimInfo* NewInfo = newCreator.createSim(outputFileName, false);
#	Line 308 | Line 369 | int main(int argc, char *argv []) {
369		Molecule* mol;
370		SimInfo::MoleculeIterator mi;
371		mol = NewInfo->beginMolecule(mi);
372	+
373		int l = 0;
374
375		for (int i = 0; i < nComponents; i++){
376		locator = new MoLocator(NewInfo->getMoleculeStamp(i),
377		NewInfo->getForceField());
378	<
379	<	if (args_info.ShellRadius_given) {
380	<	for (int n = 0; n < sites.size(); n++) {
381	<	if (componentFromSite[n] == i) {
382	<	mol = NewInfo->getMoleculeByGlobalIndex(l);
383	<	locator->placeMol(sites[n], orientations[n], mol);
384	<	l++;
378	>
379	>	if (!args_info.molFraction_given) {
380	>	for (unsigned int n = 0; n < sites.size(); n++) {
381	>	if (!isVacancy[n]) {
382	>	if (componentFromSite[n] == i) {
383	>	mol = NewInfo->getMoleculeByGlobalIndex(l);
384	>	locator->placeMol(sites[n], orientations[n], mol);
385	>	l++;
386	>	}
387		}
388		}
389		} else {
#	Line 329 | Line 393 | int main(int argc, char *argv []) {
393		l++;
394		}
395		}
396	<	}
396	>	}
397
398		//fill Hmat
399		hmat(0, 0)=  10.0*particleRadius;
#	Line 352 | Line 416 | int main(int argc, char *argv []) {
416		writer = new DumpWriter(NewInfo, outputFileName);
417
418		if (writer == NULL) {
419	<	sprintf(painCave.errMsg, "Error in creating dumpwrite object ");
419	>	sprintf(painCave.errMsg, "Error in creating dumpwriter object ");
420		painCave.isFatal = 1;
421		simError();
422		}
#	Line 364 | Line 428 | int main(int argc, char *argv []) {
428		delete writer;
429
430		// cleanup a by calling sim error.....
431	<	sprintf(painCave.errMsg, "A new OOPSE MD file called \"%s\" has been "
431	>	sprintf(painCave.errMsg, "A new OpenMD file called \"%s\" has been "
432		"generated.\n", outputFileName.c_str());
433		painCave.isFatal = 0;
434	+	painCave.severity = OPENMD_INFO;
435		simError();
436		return 0;
437		}
#	Line 382 | Line 447 | void createMdFile(const std::string&oldMdFileName,
447		//create new .md file based on old .md file
448		oldMdFile.open(oldMdFileName.c_str());
449		newMdFile.open(newMdFileName.c_str());
385	–
450		oldMdFile.getline(buffer, MAXLEN);
451	<
452	<	int i = 0;
451	>
452	>	unsigned int i = 0;
453		while (!oldMdFile.eof()) {
454	<
454	>
455		//correct molecule number
456		if (strstr(buffer, "nMol") != NULL) {
457		if(i<nMol.size()){
458	<	sprintf(buffer, "\tnMol = %i;", nMol.at(i));
458	>	sprintf(buffer, "\tnMol = %i;", nMol.at(i));
459		newMdFile << buffer << std::endl;
460		i++;
461		}
#	Line 403 | Line 467 | void createMdFile(const std::string&oldMdFileName,
467
468		oldMdFile.close();
469		newMdFile.close();
470	+
471	+	if (i != nMol.size()) {
472	+	sprintf(painCave.errMsg, "Couldn't replace the correct number of nMol\n"
473	+	"\tstatements in component blocks.  Make sure that all\n"
474	+	"\tcomponents in the template file have nMol=1");
475	+	painCave.isFatal = 1;
476	+	simError();
477	+	}
478	+
479		}
480

Comparing trunk/src/applications/nanoparticleBuilder/nanoparticleBuilder.cpp (property svn:keywords):
Revision 1075 by gezelter, Tue Oct 17 17:51:52 2006 UTC vs.
Revision 1977 by gezelter, Wed Mar 12 21:35:23 2014 UTC

#	Line 0 | Line 1
1	+	Author Id Revision Date

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines