ViewVC Help

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

root/OpenMD/branches/heatflux/src/UseTheForce/ForceField.cpp

Comparing trunk/src/UseTheForce/ForceField.cpp (file contents):
Revision 246 by gezelter, Wed Jan 12 22:41:40 2005 UTC vs.
Revision 1596 by gezelter, Mon Jul 25 17:30:53 2011 UTC

#	Line 1 | Line 1
1	<	/*
1	>	/*
2		* Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3		*
4		* The University of Notre Dame grants you ("Licensee") a
#	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, 24107 (2008).
39	+	* [4]  Vardeman & Gezelter, in progress (2009).
40		*/
41
42	<	/**
43	<	* @file ForceField.cpp
44	<	* @author tlin
45	<	* @date 11/04/2004
46	<	* @time 22:51am
47	<	* @version 1.0
48	<	*/
42	>	/**
43	>	* @file ForceField.cpp
44	>	* @author tlin
45	>	* @date 11/04/2004
46	>	* @time 22:51am
47	>	* @version 1.0
48	>	*/
49
50	+	#include <algorithm>
51		#include "UseTheForce/ForceField.hpp"
52		#include "utils/simError.h"
53	<	namespace oopse {
53	>	#include "utils/Tuple.hpp"
54	>	#include "UseTheForce/DarkSide/atype_interface.h"
55	>	#include "UseTheForce/DarkSide/fForceOptions_interface.h"
56	>	#include "UseTheForce/DarkSide/switcheroo_interface.h"
57	>	namespace OpenMD {
58
59	<	ForceField::ForceField() {
59	>	ForceField::ForceField() {
60	>
61		char* tempPath;
62		tempPath = getenv("FORCE_PARAM_PATH");
63	<
63	>
64		if (tempPath == NULL) {
65	<	//convert a macro from compiler to a string in c++
66	<	STR_DEFINE(ffPath_, FRC_PATH );
65	>	//convert a macro from compiler to a string in c++
66	>	STR_DEFINE(ffPath_, FRC_PATH );
67		} else {
68	<	ffPath_ = tempPath;
68	>	ffPath_ = tempPath;
69		}
70	<	}
70	>	}
71
72	<	AtomType* ForceField::getAtomType(const std::string &at) {
72	>
73	>	ForceField::~ForceField() {
74	>	deleteAtypes();
75	>	deleteSwitch();
76	>	}
77	>
78	>	AtomType* ForceField::getAtomType(const std::string &at) {
79		std::vector<std::string> keys;
80		keys.push_back(at);
81		return atomTypeCont_.find(keys);
82	<	}
82	>	}
83
84	<	BondType* ForceField::getBondType(const std::string &at1, const std::string &at2) {
84	>	BondType* ForceField::getBondType(const std::string &at1,
85	>	const std::string &at2) {
86		std::vector<std::string> keys;
87		keys.push_back(at1);
88		keys.push_back(at2);
#	Line 77 | Line 90 | BondType* ForceField::getBondType(const std::string &a
90		//try exact match first
91		BondType* bondType = bondTypeCont_.find(keys);
92		if (bondType) {
93	<	return bondType;
93	>	return bondType;
94		} else {
95	<	//if no exact match found, try wild card match
96	<	return bondTypeCont_.find(keys, wildCardAtomTypeName_);
97	<	}
95	>	AtomType* atype1;
96	>	AtomType* atype2;
97	>	std::vector<std::string> at1key;
98	>	at1key.push_back(at1);
99	>	atype1 = atomTypeCont_.find(at1key);
100	>
101	>	std::vector<std::string> at2key;
102	>	at2key.push_back(at2);
103	>	atype2 = atomTypeCont_.find(at2key);
104
105	<	}
105	>	// query atom types for their chains of responsibility
106	>	std::vector<AtomType*> at1Chain = atype1->allYourBase();
107	>	std::vector<AtomType*> at2Chain = atype2->allYourBase();
108
109	<	BendType* ForceField::getBendType(const std::string &at1, const std::string &at2,
110	<	const std::string &at3) {
109	>	std::vector<AtomType*>::iterator i;
110	>	std::vector<AtomType*>::iterator j;
111	>
112	>	int ii = 0;
113	>	int jj = 0;
114	>	int bondTypeScore;
115	>
116	>	std::vector<std::pair<int, std::vector<std::string> > > foundBonds;
117	>
118	>	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
119	>	jj = 0;
120	>	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
121	>
122	>	bondTypeScore = ii + jj;
123	>
124	>	std::vector<std::string> myKeys;
125	>	myKeys.push_back((*i)->getName());
126	>	myKeys.push_back((*j)->getName());
127	>
128	>	BondType* bondType = bondTypeCont_.find(myKeys);
129	>	if (bondType) {
130	>	foundBonds.push_back(std::make_pair(bondTypeScore, myKeys));
131	>	}
132	>	jj++;
133	>	}
134	>	ii++;
135	>	}
136	>
137	>
138	>	if (foundBonds.size() > 0) {
139	>	// sort the foundBonds by the score:
140	>	std::sort(foundBonds.begin(), foundBonds.end());
141	>
142	>	int bestScore = foundBonds[0].first;
143	>	std::vector<std::string> theKeys = foundBonds[0].second;
144	>
145	>	BondType* bestType = bondTypeCont_.find(theKeys);
146	>
147	>	return bestType;
148	>	} else {
149	>	//if no exact match found, try wild card match
150	>	return bondTypeCont_.find(keys, wildCardAtomTypeName_);
151	>	}
152	>	}
153	>	}
154	>
155	>	BendType* ForceField::getBendType(const std::string &at1,
156	>	const std::string &at2,
157	>	const std::string &at3) {
158		std::vector<std::string> keys;
159		keys.push_back(at1);
160		keys.push_back(at2);
#	Line 95 | Line 163 | BendType* ForceField::getBendType(const std::string &a
163		//try exact match first
164		BendType* bendType = bendTypeCont_.find(keys);
165		if (bendType) {
166	<	return bendType;
166	>	return bendType;
167		} else {
100	–	//if no exact match found, try wild card match
101	–	return bendTypeCont_.find(keys, wildCardAtomTypeName_);
102	–	}
103	–	}
168
169	<	TorsionType* ForceField::getTorsionType(const std::string &at1, const std::string &at2,
170	<	const std::string &at3, const std::string &at4) {
169	>	AtomType* atype1;
170	>	AtomType* atype2;
171	>	AtomType* atype3;
172	>	std::vector<std::string> at1key;
173	>	at1key.push_back(at1);
174	>	atype1 = atomTypeCont_.find(at1key);
175	>
176	>	std::vector<std::string> at2key;
177	>	at2key.push_back(at2);
178	>	atype2 = atomTypeCont_.find(at2key);
179	>
180	>	std::vector<std::string> at3key;
181	>	at3key.push_back(at3);
182	>	atype3 = atomTypeCont_.find(at3key);
183	>
184	>	// query atom types for their chains of responsibility
185	>	std::vector<AtomType*> at1Chain = atype1->allYourBase();
186	>	std::vector<AtomType*> at2Chain = atype2->allYourBase();
187	>	std::vector<AtomType*> at3Chain = atype3->allYourBase();
188	>
189	>	std::vector<AtomType*>::iterator i;
190	>	std::vector<AtomType*>::iterator j;
191	>	std::vector<AtomType*>::iterator k;
192	>
193	>	int ii = 0;
194	>	int jj = 0;
195	>	int kk = 0;
196	>	int IKscore;
197	>
198	>	std::vector<tuple3<int, int, std::vector<std::string> > > foundBends;
199	>
200	>	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
201	>	ii = 0;
202	>	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
203	>	kk = 0;
204	>	for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
205	>
206	>	IKscore = ii + kk;
207	>
208	>	std::vector<std::string> myKeys;
209	>	myKeys.push_back((*i)->getName());
210	>	myKeys.push_back((*j)->getName());
211	>	myKeys.push_back((*k)->getName());
212	>
213	>	BendType* bendType = bendTypeCont_.find(myKeys);
214	>	if (bendType) {
215	>	foundBends.push_back( make_tuple3(jj, IKscore, myKeys) );
216	>	}
217	>	kk++;
218	>	}
219	>	ii++;
220	>	}
221	>	jj++;
222	>	}
223	>
224	>	if (foundBends.size() > 0) {
225	>	std::sort(foundBends.begin(), foundBends.end());
226	>	int jscore = foundBends[0].first;
227	>	int ikscore = foundBends[0].second;
228	>	std::vector<std::string> theKeys = foundBends[0].third;
229	>
230	>	BendType* bestType = bendTypeCont_.find(theKeys);
231	>	return bestType;
232	>	} else {
233	>	//if no exact match found, try wild card match
234	>	return bendTypeCont_.find(keys, wildCardAtomTypeName_);
235	>	}
236	>	}
237	>	}
238	>
239	>	TorsionType* ForceField::getTorsionType(const std::string &at1,
240	>	const std::string &at2,
241	>	const std::string &at3,
242	>	const std::string &at4) {
243		std::vector<std::string> keys;
244		keys.push_back(at1);
245		keys.push_back(at2);
246		keys.push_back(at3);
247		keys.push_back(at4);
248
249	+
250	+	//try exact match first
251		TorsionType* torsionType = torsionTypeCont_.find(keys);
252		if (torsionType) {
253	<	return torsionType;
253	>	return torsionType;
254		} else {
117	–	//if no exact match found, try wild card match
118	–	return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
119	–	}
120	–
121	–	return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
255
256	<	}
256	>	AtomType* atype1;
257	>	AtomType* atype2;
258	>	AtomType* atype3;
259	>	AtomType* atype4;
260	>	std::vector<std::string> at1key;
261	>	at1key.push_back(at1);
262	>	atype1 = atomTypeCont_.find(at1key);
263	>
264	>	std::vector<std::string> at2key;
265	>	at2key.push_back(at2);
266	>	atype2 = atomTypeCont_.find(at2key);
267
268	<	BondType* ForceField::getExactBondType(const std::string &at1, const std::string &at2){
268	>	std::vector<std::string> at3key;
269	>	at3key.push_back(at3);
270	>	atype3 = atomTypeCont_.find(at3key);
271	>
272	>	std::vector<std::string> at4key;
273	>	at4key.push_back(at4);
274	>	atype4 = atomTypeCont_.find(at4key);
275	>
276	>	// query atom types for their chains of responsibility
277	>	std::vector<AtomType*> at1Chain = atype1->allYourBase();
278	>	std::vector<AtomType*> at2Chain = atype2->allYourBase();
279	>	std::vector<AtomType*> at3Chain = atype3->allYourBase();
280	>	std::vector<AtomType*> at4Chain = atype4->allYourBase();
281	>
282	>	std::vector<AtomType*>::iterator i;
283	>	std::vector<AtomType*>::iterator j;
284	>	std::vector<AtomType*>::iterator k;
285	>	std::vector<AtomType*>::iterator l;
286	>
287	>	int ii = 0;
288	>	int jj = 0;
289	>	int kk = 0;
290	>	int ll = 0;
291	>	int ILscore;
292	>	int JKscore;
293	>
294	>	std::vector<tuple3<int, int, std::vector<std::string> > > foundTorsions;
295	>
296	>	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
297	>	kk = 0;
298	>	for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
299	>	ii = 0;
300	>	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
301	>	ll = 0;
302	>	for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
303	>
304	>	ILscore = ii + ll;
305	>	JKscore = jj + kk;
306	>
307	>	std::vector<std::string> myKeys;
308	>	myKeys.push_back((*i)->getName());
309	>	myKeys.push_back((*j)->getName());
310	>	myKeys.push_back((*k)->getName());
311	>	myKeys.push_back((*l)->getName());
312	>
313	>	TorsionType* torsionType = torsionTypeCont_.find(myKeys);
314	>	if (torsionType) {
315	>	foundTorsions.push_back( make_tuple3(JKscore, ILscore, myKeys) );
316	>	}
317	>	ll++;
318	>	}
319	>	ii++;
320	>	}
321	>	kk++;
322	>	}
323	>	jj++;
324	>	}
325	>
326	>	if (foundTorsions.size() > 0) {
327	>	std::sort(foundTorsions.begin(), foundTorsions.end());
328	>	int jkscore = foundTorsions[0].first;
329	>	int ilscore = foundTorsions[0].second;
330	>	std::vector<std::string> theKeys = foundTorsions[0].third;
331	>
332	>	TorsionType* bestType = torsionTypeCont_.find(theKeys);
333	>	return bestType;
334	>	} else {
335	>	//if no exact match found, try wild card match
336	>	return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
337	>	}
338	>	}
339	>	}
340	>
341	>	InversionType* ForceField::getInversionType(const std::string &at1,
342	>	const std::string &at2,
343	>	const std::string &at3,
344	>	const std::string &at4) {
345		std::vector<std::string> keys;
346		keys.push_back(at1);
347		keys.push_back(at2);
348	<	return bondTypeCont_.find(keys);
349	<	}
348	>	keys.push_back(at3);
349	>	keys.push_back(at4);
350
351	<	BendType* ForceField::getExactBendType(const std::string &at1, const std::string &at2,
352	<	const std::string &at3){
351	>	//try exact match first
352	>	InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(keys);
353	>	if (inversionType) {
354	>	return inversionType;
355	>	} else {
356	>
357	>	AtomType* atype1;
358	>	AtomType* atype2;
359	>	AtomType* atype3;
360	>	AtomType* atype4;
361	>	std::vector<std::string> at1key;
362	>	at1key.push_back(at1);
363	>	atype1 = atomTypeCont_.find(at1key);
364	>
365	>	std::vector<std::string> at2key;
366	>	at2key.push_back(at2);
367	>	atype2 = atomTypeCont_.find(at2key);
368	>
369	>	std::vector<std::string> at3key;
370	>	at3key.push_back(at3);
371	>	atype3 = atomTypeCont_.find(at3key);
372	>
373	>	std::vector<std::string> at4key;
374	>	at4key.push_back(at4);
375	>	atype4 = atomTypeCont_.find(at4key);
376	>
377	>	// query atom types for their chains of responsibility
378	>	std::vector<AtomType*> at1Chain = atype1->allYourBase();
379	>	std::vector<AtomType*> at2Chain = atype2->allYourBase();
380	>	std::vector<AtomType*> at3Chain = atype3->allYourBase();
381	>	std::vector<AtomType*> at4Chain = atype4->allYourBase();
382	>
383	>	std::vector<AtomType*>::iterator i;
384	>	std::vector<AtomType*>::iterator j;
385	>	std::vector<AtomType*>::iterator k;
386	>	std::vector<AtomType*>::iterator l;
387	>
388	>	int ii = 0;
389	>	int jj = 0;
390	>	int kk = 0;
391	>	int ll = 0;
392	>	int Iscore;
393	>	int JKLscore;
394	>
395	>	std::vector<tuple3<int, int, std::vector<std::string> > > foundInversions;
396	>
397	>	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
398	>	kk = 0;
399	>	for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
400	>	ii = 0;
401	>	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
402	>	ll = 0;
403	>	for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
404	>
405	>	Iscore = ii;
406	>	JKLscore = jj + kk + ll;
407	>
408	>	std::vector<std::string> myKeys;
409	>	myKeys.push_back((*i)->getName());
410	>	myKeys.push_back((*j)->getName());
411	>	myKeys.push_back((*k)->getName());
412	>	myKeys.push_back((*l)->getName());
413	>
414	>	InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(myKeys);
415	>	if (inversionType) {
416	>	foundInversions.push_back( make_tuple3(Iscore, JKLscore, myKeys) );
417	>	}
418	>	ll++;
419	>	}
420	>	ii++;
421	>	}
422	>	kk++;
423	>	}
424	>	jj++;
425	>	}
426	>
427	>	if (foundInversions.size() > 0) {
428	>	std::sort(foundInversions.begin(), foundInversions.end());
429	>	int iscore = foundInversions[0].first;
430	>	int jklscore = foundInversions[0].second;
431	>	std::vector<std::string> theKeys = foundInversions[0].third;
432	>
433	>	InversionType* bestType = inversionTypeCont_.permutedFindSkippingFirstElement(theKeys);
434	>	return bestType;
435	>	} else {
436	>	//if no exact match found, try wild card match
437	>	return inversionTypeCont_.find(keys, wildCardAtomTypeName_);
438	>	}
439	>	}
440	>	}
441	>
442	>	NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
443		std::vector<std::string> keys;
444		keys.push_back(at1);
445		keys.push_back(at2);
446	+
447	+	//try exact match first
448	+	NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(keys);
449	+	if (nbiType) {
450	+	return nbiType;
451	+	} else {
452	+	AtomType* atype1;
453	+	AtomType* atype2;
454	+	std::vector<std::string> at1key;
455	+	at1key.push_back(at1);
456	+	atype1 = atomTypeCont_.find(at1key);
457	+
458	+	std::vector<std::string> at2key;
459	+	at2key.push_back(at2);
460	+	atype2 = atomTypeCont_.find(at2key);
461	+
462	+	// query atom types for their chains of responsibility
463	+	std::vector<AtomType*> at1Chain = atype1->allYourBase();
464	+	std::vector<AtomType*> at2Chain = atype2->allYourBase();
465	+
466	+	std::vector<AtomType*>::iterator i;
467	+	std::vector<AtomType*>::iterator j;
468	+
469	+	int ii = 0;
470	+	int jj = 0;
471	+	int nbiTypeScore;
472	+
473	+	std::vector<std::pair<int, std::vector<std::string> > > foundNBI;
474	+
475	+	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
476	+	jj = 0;
477	+	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
478	+
479	+	nbiTypeScore = ii + jj;
480	+
481	+	std::vector<std::string> myKeys;
482	+	myKeys.push_back((*i)->getName());
483	+	myKeys.push_back((*j)->getName());
484	+
485	+	NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(myKeys);
486	+	if (nbiType) {
487	+	foundNBI.push_back(std::make_pair(nbiTypeScore, myKeys));
488	+	}
489	+	jj++;
490	+	}
491	+	ii++;
492	+	}
493	+
494	+
495	+	if (foundNBI.size() > 0) {
496	+	// sort the foundNBI by the score:
497	+	std::sort(foundNBI.begin(), foundNBI.end());
498	+
499	+	int bestScore = foundNBI[0].first;
500	+	std::vector<std::string> theKeys = foundNBI[0].second;
501	+
502	+	NonBondedInteractionType* bestType = nonBondedInteractionTypeCont_.find(theKeys);
503	+	return bestType;
504	+	} else {
505	+	//if no exact match found, try wild card match
506	+	return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
507	+	}
508	+	}
509	+	}
510	+
511	+	BondType* ForceField::getExactBondType(const std::string &at1,
512	+	const std::string &at2){
513	+	std::vector<std::string> keys;
514	+	keys.push_back(at1);
515	+	keys.push_back(at2);
516	+	return bondTypeCont_.find(keys);
517	+	}
518	+
519	+	BendType* ForceField::getExactBendType(const std::string &at1,
520	+	const std::string &at2,
521	+	const std::string &at3){
522	+	std::vector<std::string> keys;
523	+	keys.push_back(at1);
524	+	keys.push_back(at2);
525		keys.push_back(at3);
526		return bendTypeCont_.find(keys);
527	<	}
528	<
529	<	TorsionType* ForceField::getExactTorsionType(const std::string &at1, const std::string &at2,
530	<	const std::string &at3, const std::string &at4){
527	>	}
528	>
529	>	TorsionType* ForceField::getExactTorsionType(const std::string &at1,
530	>	const std::string &at2,
531	>	const std::string &at3,
532	>	const std::string &at4){
533		std::vector<std::string> keys;
534		keys.push_back(at1);
535		keys.push_back(at2);
536		keys.push_back(at3);
537		keys.push_back(at4);
538		return torsionTypeCont_.find(keys);
539	<	}
540	<	bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
539	>	}
540	>
541	>	InversionType* ForceField::getExactInversionType(const std::string &at1,
542	>	const std::string &at2,
543	>	const std::string &at3,
544	>	const std::string &at4){
545		std::vector<std::string> keys;
546	+	keys.push_back(at1);
547	+	keys.push_back(at2);
548	+	keys.push_back(at3);
549	+	keys.push_back(at4);
550	+	return inversionTypeCont_.find(keys);
551	+	}
552	+
553	+	NonBondedInteractionType* ForceField::getExactNonBondedInteractionType(const std::string &at1, const std::string &at2){
554	+	std::vector<std::string> keys;
555	+	keys.push_back(at1);
556	+	keys.push_back(at2);
557	+	return nonBondedInteractionTypeCont_.find(keys);
558	+	}
559	+
560	+
561	+	bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
562	+	std::vector<std::string> keys;
563		keys.push_back(at);
564		return atomTypeCont_.add(keys, atomType);
565	<	}
565	>	}
566
567	<	bool ForceField::addBondType(const std::string &at1, const std::string &at2, BondType* bondType) {
567	>	bool ForceField::replaceAtomType(const std::string &at, AtomType* atomType) {
568		std::vector<std::string> keys;
569	+	keys.push_back(at);
570	+	return atomTypeCont_.replace(keys, atomType);
571	+	}
572	+
573	+	bool ForceField::addBondType(const std::string &at1, const std::string &at2,
574	+	BondType* bondType) {
575	+	std::vector<std::string> keys;
576		keys.push_back(at1);
577		keys.push_back(at2);
578	<	return bondTypeCont_.add(keys, bondType);
579	<
580	<	}
581	<
582	<	bool ForceField::addBendType(const std::string &at1, const std::string &at2,
165	<	const std::string &at3, BendType* bendType) {
578	>	return bondTypeCont_.add(keys, bondType);
579	>	}
580	>
581	>	bool ForceField::addBendType(const std::string &at1, const std::string &at2,
582	>	const std::string &at3, BendType* bendType) {
583		std::vector<std::string> keys;
584		keys.push_back(at1);
585		keys.push_back(at2);
586		keys.push_back(at3);
587		return bendTypeCont_.add(keys, bendType);
588	<	}
589	<
590	<	bool ForceField::addTorsionType(const std::string &at1, const std::string &at2,
591	<	const std::string &at3, const std::string &at4, TorsionType* torsionType) {
588	>	}
589	>
590	>	bool ForceField::addTorsionType(const std::string &at1,
591	>	const std::string &at2,
592	>	const std::string &at3,
593	>	const std::string &at4,
594	>	TorsionType* torsionType) {
595		std::vector<std::string> keys;
596		keys.push_back(at1);
597		keys.push_back(at2);
598		keys.push_back(at3);
599		keys.push_back(at4);
600		return torsionTypeCont_.add(keys, torsionType);
601	<	}
601	>	}
602
603	<	double ForceField::getRcutFromAtomType(AtomType* at) {
603	>	bool ForceField::addInversionType(const std::string &at1,
604	>	const std::string &at2,
605	>	const std::string &at3,
606	>	const std::string &at4,
607	>	InversionType* inversionType) {
608	>	std::vector<std::string> keys;
609	>	keys.push_back(at1);
610	>	keys.push_back(at2);
611	>	keys.push_back(at3);
612	>	keys.push_back(at4);
613	>	return inversionTypeCont_.add(keys, inversionType);
614	>	}
615	>
616	>	bool ForceField::addNonBondedInteractionType(const std::string &at1,
617	>	const std::string &at2,
618	>	NonBondedInteractionType* nbiType) {
619	>	std::vector<std::string> keys;
620	>	keys.push_back(at1);
621	>	keys.push_back(at2);
622	>	return nonBondedInteractionTypeCont_.add(keys, nbiType);
623	>	}
624	>
625	>	RealType ForceField::getRcutFromAtomType(AtomType* at) {
626		/**@todo */
627		GenericData* data;
628	<	double rcut = 0.0;
629	<
628	>	RealType rcut = 0.0;
629	>
630		if (at->isLennardJones()) {
631	<	data = at->getPropertyByName("LennardJones");
632	<	if (data != NULL) {
633	<	LJParamGenericData* ljData = dynamic_cast<LJParamGenericData*>(data);
634	<
635	<	if (ljData != NULL) {
636	<	LJParam ljParam = ljData->getData();
637	<
638	<	//by default use 2.5*sigma as cutoff radius
639	<	rcut = 2.5 * ljParam.sigma;
640	<
641	<	} else {
642	<	sprintf( painCave.errMsg,
643	<	"Can not cast GenericData to LJParam\n");
644	<	painCave.severity = OOPSE_ERROR;
645	<	painCave.isFatal = 1;
646	<	simError();
647	<	}
648	<	} else {
649	<	sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
650	<	painCave.severity = OOPSE_ERROR;
651	<	painCave.isFatal = 1;
652	<	simError();
653	<	}
631	>	data = at->getPropertyByName("LennardJones");
632	>	if (data != NULL) {
633	>	LJParamGenericData* ljData = dynamic_cast<LJParamGenericData*>(data);
634	>
635	>	if (ljData != NULL) {
636	>	LJParam ljParam = ljData->getData();
637	>
638	>	//by default use 2.5*sigma as cutoff radius
639	>	rcut = 2.5 * ljParam.sigma;
640	>
641	>	} else {
642	>	sprintf( painCave.errMsg,
643	>	"Can not cast GenericData to LJParam\n");
644	>	painCave.severity = OPENMD_ERROR;
645	>	painCave.isFatal = 1;
646	>	simError();
647	>	}
648	>	} else {
649	>	sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
650	>	painCave.severity = OPENMD_ERROR;
651	>	painCave.isFatal = 1;
652	>	simError();
653	>	}
654		}
213	–
655		return rcut;
656	<	}
656	>	}
657	>
658
659	<
218	<	ifstrstream* ForceField::openForceFieldFile(const std::string& filename) {
659	>	ifstrstream* ForceField::openForceFieldFile(const std::string& filename) {
660		std::string forceFieldFilename(filename);
661		ifstrstream* ffStream = new ifstrstream();
662
#	Line 223 | Line 664 | ifstrstream* ForceField::openForceFieldFile(const std:
664		ffStream->open(forceFieldFilename.c_str());
665		if(!ffStream->is_open()){
666
667	<	forceFieldFilename = ffPath_ + "/" + forceFieldFilename;
668	<	ffStream->open( forceFieldFilename.c_str() );
667	>	forceFieldFilename = ffPath_ + "/" + forceFieldFilename;
668	>	ffStream->open( forceFieldFilename.c_str() );
669
670	<	//if current directory does not contain the force field file,
671	<	//try to open it in the path
672	<	if(!ffStream->is_open()){
670	>	//if current directory does not contain the force field file,
671	>	//try to open it in the path
672	>	if(!ffStream->is_open()){
673
674	<	sprintf( painCave.errMsg,
675	<	"Error opening the force field parameter file:\n"
676	<	"\t%s\n"
677	<	"\tHave you tried setting the FORCE_PARAM_PATH environment "
678	<	"variable?\n",
679	<	forceFieldFilename.c_str() );
680	<	painCave.severity = OOPSE_ERROR;
681	<	painCave.isFatal = 1;
682	<	simError();
683	<	}
674	>	sprintf( painCave.errMsg,
675	>	"Error opening the force field parameter file:\n"
676	>	"\t%s\n"
677	>	"\tHave you tried setting the FORCE_PARAM_PATH environment "
678	>	"variable?\n",
679	>	forceFieldFilename.c_str() );
680	>	painCave.severity = OPENMD_ERROR;
681	>	painCave.isFatal = 1;
682	>	simError();
683	>	}
684		}
244	–
685		return ffStream;
686	+	}
687
688	<	}
689	<
690	<	} //end namespace oopse
688	>	void ForceField::setFortranForceOptions(){
689	>	ForceOptions theseFortranOptions;
690	>	forceFieldOptions_.makeFortranOptions(theseFortranOptions);
691	>	setfForceOptions(&theseFortranOptions);
692	>	}
693	>	} //end namespace OpenMD

Comparing trunk/src/UseTheForce/ForceField.cpp (property svn:keywords):
Revision 246 by gezelter, Wed Jan 12 22:41:40 2005 UTC vs.
Revision 1596 by gezelter, Mon Jul 25 17:30:53 2011 UTC

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

Diff Legend

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