ViewVC Help

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

root/OpenMD/trunk/src/brains/ForceField.cpp

Comparing:
trunk/src/UseTheForce/ForceField.cpp (file contents), Revision 1269 by gezelter, Tue Jul 1 13:28:23 2008 UTC vs.
trunk/src/brains/ForceField.cpp (file contents), Revision 1880 by gezelter, Mon Jun 17 18:28:30 2013 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		/**
44		* @file ForceField.cpp
45		* @author tlin
46		* @date 11/04/2004
46	–	* @time 22:51am
47		* @version 1.0
48		*/
49
50		#include <algorithm>
51	<	#include "UseTheForce/ForceField.hpp"
51	>	#include "brains/ForceField.hpp"
52		#include "utils/simError.h"
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 oopse {
53
54	<	ForceField::ForceField() {
54	>	#include "io/OptionSectionParser.hpp"
55	>	#include "io/BaseAtomTypesSectionParser.hpp"
56	>	#include "io/DirectionalAtomTypesSectionParser.hpp"
57	>	#include "io/AtomTypesSectionParser.hpp"
58	>	#include "io/BendTypesSectionParser.hpp"
59	>	#include "io/BondTypesSectionParser.hpp"
60	>	#include "io/ChargeAtomTypesSectionParser.hpp"
61	>	#include "io/EAMAtomTypesSectionParser.hpp"
62	>	#include "io/FluctuatingChargeAtomTypesSectionParser.hpp"
63	>	#include "io/GayBerneAtomTypesSectionParser.hpp"
64	>	#include "io/InversionTypesSectionParser.hpp"
65	>	#include "io/LennardJonesAtomTypesSectionParser.hpp"
66	>	#include "io/MultipoleAtomTypesSectionParser.hpp"
67	>	#include "io/NonBondedInteractionsSectionParser.hpp"
68	>	#include "io/PolarizableAtomTypesSectionParser.hpp"
69	>	#include "io/SCAtomTypesSectionParser.hpp"
70	>	#include "io/ShapeAtomTypesSectionParser.hpp"
71	>	#include "io/StickyAtomTypesSectionParser.hpp"
72	>	#include "io/StickyPowerAtomTypesSectionParser.hpp"
73	>	#include "io/TorsionTypesSectionParser.hpp"
74	>
75	>	#include "types/LennardJonesAdapter.hpp"
76	>	#include "types/EAMAdapter.hpp"
77	>	#include "types/SuttonChenAdapter.hpp"
78	>	#include "types/GayBerneAdapter.hpp"
79	>	#include "types/StickyAdapter.hpp"
80	>
81	>	namespace OpenMD {
82	>
83	>	ForceField::ForceField(std::string ffName) {
84	>
85		char* tempPath;
86		tempPath = getenv("FORCE_PARAM_PATH");
87	<
87	>
88		if (tempPath == NULL) {
89		//convert a macro from compiler to a string in c++
90		STR_DEFINE(ffPath_, FRC_PATH );
91		} else {
92		ffPath_ = tempPath;
93		}
94	+
95	+	setForceFieldFileName(ffName + ".frc");
96	+
97	+	/**
98	+	* The order of adding section parsers is important.
99	+	*
100	+	* OptionSectionParser must come first to set options for other
101	+	* parsers
102	+	*
103	+	* DirectionalAtomTypesSectionParser should be added before
104	+	* AtomTypesSectionParser, and these two section parsers will
105	+	* actually create "real" AtomTypes (AtomTypesSectionParser will
106	+	* create AtomType and DirectionalAtomTypesSectionParser will
107	+	* create DirectionalAtomType, which is a subclass of AtomType and
108	+	* should come first).
109	+	*
110	+	* Other AtomTypes Section Parsers will not create the "real"
111	+	* AtomType, they only add and set some attributes of the AtomType
112	+	* (via the Adapters). Thus ordering of these is not important.
113	+	* AtomTypesSectionParser should be added before other atom type
114	+	*
115	+	* The order of BondTypesSectionParser, BendTypesSectionParser and
116	+	* TorsionTypesSectionParser, etc. are not important.
117	+	*/
118	+
119	+	spMan_.push_back(new OptionSectionParser(forceFieldOptions_));
120	+	spMan_.push_back(new BaseAtomTypesSectionParser());
121	+	spMan_.push_back(new DirectionalAtomTypesSectionParser(forceFieldOptions_));
122	+	spMan_.push_back(new AtomTypesSectionParser());
123	+
124	+	spMan_.push_back(new LennardJonesAtomTypesSectionParser(forceFieldOptions_));
125	+	spMan_.push_back(new ChargeAtomTypesSectionParser(forceFieldOptions_));
126	+	spMan_.push_back(new MultipoleAtomTypesSectionParser(forceFieldOptions_));
127	+	spMan_.push_back(new FluctuatingChargeAtomTypesSectionParser(forceFieldOptions_));
128	+	spMan_.push_back(new PolarizableAtomTypesSectionParser(forceFieldOptions_));
129	+	spMan_.push_back(new GayBerneAtomTypesSectionParser(forceFieldOptions_));
130	+	spMan_.push_back(new EAMAtomTypesSectionParser(forceFieldOptions_));
131	+	spMan_.push_back(new SCAtomTypesSectionParser(forceFieldOptions_));
132	+	spMan_.push_back(new ShapeAtomTypesSectionParser(forceFieldOptions_));
133	+	spMan_.push_back(new StickyAtomTypesSectionParser(forceFieldOptions_));
134	+	spMan_.push_back(new StickyPowerAtomTypesSectionParser(forceFieldOptions_));
135	+
136	+	spMan_.push_back(new BondTypesSectionParser(forceFieldOptions_));
137	+	spMan_.push_back(new BendTypesSectionParser(forceFieldOptions_));
138	+	spMan_.push_back(new TorsionTypesSectionParser(forceFieldOptions_));
139	+	spMan_.push_back(new InversionTypesSectionParser(forceFieldOptions_));
140	+
141	+	spMan_.push_back(new NonBondedInteractionsSectionParser(forceFieldOptions_));
142		}
143
144	+	void ForceField::parse(const std::string& filename) {
145	+	ifstrstream* ffStream;
146
147	<	ForceField::~ForceField() {
148	<	deleteAtypes();
149	<	deleteSwitch();
147	>	ffStream = openForceFieldFile(filename);
148	>
149	>	spMan_.parse(*ffStream, *this);
150	>
151	>	ForceField::AtomTypeContainer::MapTypeIterator i;
152	>	AtomType* at;
153	>
154	>	for (at = atomTypeCont_.beginType(i); at != NULL;
155	>	at = atomTypeCont_.nextType(i)) {
156	>
157	>	// useBase sets the responsibilities, and these have to be done
158	>	// after the atomTypes and Base types have all been scanned:
159	>
160	>	std::vector<AtomType*> ayb = at->allYourBase();
161	>	if (ayb.size() > 1) {
162	>	for (int j = ayb.size()-1; j > 0; j--) {
163	>
164	>	ayb[j-1]->useBase(ayb[j]);
165	>
166	>	}
167	>	}
168	>	}
169	>
170	>	delete ffStream;
171		}
172
173	+	/**
174	+	* getAtomType by string
175	+	*
176	+	* finds the requested atom type in this force field using the string
177	+	* name of the atom type.
178	+	*/
179		AtomType* ForceField::getAtomType(const std::string &at) {
180		std::vector<std::string> keys;
181		keys.push_back(at);
182		return atomTypeCont_.find(keys);
183		}
184
185	+	/**
186	+	* getAtomType by ident
187	+	*
188	+	* finds the requested atom type in this force field using the
189	+	* integer ident instead of the string name of the atom type.
190	+	*/
191	+	AtomType* ForceField::getAtomType(int ident) {
192	+	std::string at = atypeIdentToName.find(ident)->second;
193	+	return getAtomType(at);
194	+	}
195	+
196		BondType* ForceField::getBondType(const std::string &at1,
197		const std::string &at2) {
198		std::vector<std::string> keys;
#	Line 114 | Line 227 | namespace oopse {
227
228		std::vector<std::pair<int, std::vector<std::string> > > foundBonds;
229
230	<	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
230	>	for (i = at1Chain.begin(); i != at1Chain.end(); ++i) {
231		jj = 0;
232	<	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
232	>	for (j = at2Chain.begin(); j != at2Chain.end(); ++j) {
233
234		bondTypeScore = ii + jj;
235
#	Line 133 | Line 246 | namespace oopse {
246		ii++;
247		}
248
136	–	// sort the foundBonds by the score:
249
250	<	std::sort(foundBonds.begin(), foundBonds.end());
251	<
252	<	int bestScore = foundBonds[0].first;
253	<	std::vector<std::string> theKeys = foundBonds[0].second;
254	<
255	<	std::cout << "best matching bond = " << theKeys[0] << "\t" << theKeys[1]  << "\t(score = "<< bestScore << ")\n";
256	<	BondType* bestType = bondTypeCont_.find(theKeys);
257	<	if (bestType)
258	<	return bestType;
259	<	else {
260	<	//if no exact match found, try wild card match
261	<	return bondTypeCont_.find(keys, wildCardAtomTypeName_);
250	>	if (!foundBonds.empty()) {
251	>	// sort the foundBonds by the score:
252	>	std::sort(foundBonds.begin(), foundBonds.end());
253	>
254	>	std::vector<std::string> theKeys = foundBonds[0].second;
255	>
256	>	BondType* bestType = bondTypeCont_.find(theKeys);
257	>
258	>	return bestType;
259	>	} else {
260	>	//if no exact match found, try wild card match
261	>	return bondTypeCont_.find(keys, wildCardAtomTypeName_);
262		}
263		}
264		}
#	Line 196 | Line 308 | namespace oopse {
308
309		std::vector<tuple3<int, int, std::vector<std::string> > > foundBends;
310
311	<	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
311	>	for (j = at2Chain.begin(); j != at2Chain.end(); ++j) {
312		ii = 0;
313	<	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
313	>	for (i = at1Chain.begin(); i != at1Chain.end(); ++i) {
314		kk = 0;
315	<	for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
315	>	for (k = at3Chain.begin(); k != at3Chain.end(); ++k) {
316
317		IKscore = ii + kk;
318
#	Line 220 | Line 332 | namespace oopse {
332		jj++;
333		}
334
335	<	std::sort(foundBends.begin(), foundBends.end());
336	<
337	<	int jscore = foundBends[0].first;
338	<	int ikscore = foundBends[0].second;
339	<	std::vector<std::string> theKeys = foundBends[0].third;
340	<
341	<	std::cout << "best matching bend = " << theKeys[0] << "\t" <<theKeys[1]  << "\t" << theKeys[2] << "\t(scores = "<< jscore << "\t" << ikscore << ")\n";
230	<
231	<	BendType* bestType = bendTypeCont_.find(theKeys);
232	<	if (bestType)
233	<	return bestType;
234	<	else {
235	<
335	>	if (!foundBends.empty()) {
336	>	std::sort(foundBends.begin(), foundBends.end());
337	>	std::vector<std::string> theKeys = foundBends[0].third;
338	>
339	>	BendType* bestType = bendTypeCont_.find(theKeys);
340	>	return bestType;
341	>	} else {
342		//if no exact match found, try wild card match
343		return bendTypeCont_.find(keys, wildCardAtomTypeName_);
344		}
345		}
346		}
347
242	–
348		TorsionType* ForceField::getTorsionType(const std::string &at1,
349		const std::string &at2,
350		const std::string &at3,
#	Line 297 | Line 402 | namespace oopse {
402
403		std::vector<tuple3<int, int, std::vector<std::string> > > foundTorsions;
404
405	<	for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
405	>	for (j = at2Chain.begin(); j != at2Chain.end(); ++j) {
406		kk = 0;
407	<	for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
407	>	for (k = at3Chain.begin(); k != at3Chain.end(); ++k) {
408		ii = 0;
409	<	for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
409	>	for (i = at1Chain.begin(); i != at1Chain.end(); ++i) {
410		ll = 0;
411	<	for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
411	>	for (l = at4Chain.begin(); l != at4Chain.end(); ++l) {
412
413		ILscore = ii + ll;
414		JKscore = jj + kk;
#	Line 327 | Line 432 | namespace oopse {
432		jj++;
433		}
434
435	<	std::sort(foundTorsions.begin(), foundTorsions.end());
435	>	if (!foundTorsions.empty()) {
436	>	std::sort(foundTorsions.begin(), foundTorsions.end());
437	>	std::vector<std::string> theKeys = foundTorsions[0].third;
438	>
439	>	TorsionType* bestType = torsionTypeCont_.find(theKeys);
440	>	return bestType;
441	>	} else {
442	>	//if no exact match found, try wild card match
443	>	return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
444	>	}
445	>	}
446	>	}
447
448	<	int jkscore = foundTorsions[0].first;
449	<	int ilscore = foundTorsions[0].second;
450	<	std::vector<std::string> theKeys = foundTorsions[0].third;
448	>	InversionType* ForceField::getInversionType(const std::string &at1,
449	>	const std::string &at2,
450	>	const std::string &at3,
451	>	const std::string &at4) {
452	>	std::vector<std::string> keys;
453	>	keys.push_back(at1);
454	>	keys.push_back(at2);
455	>	keys.push_back(at3);
456	>	keys.push_back(at4);
457
458	<	std::cout << "best matching torsion = " << theKeys[0] << "\t" <<theKeys[1]  << "\t" << theKeys[2] << "\t" << theKeys[3] << "\t(scores = "<< jkscore << "\t" << ilscore << ")\n";
458	>	//try exact match first
459	>	InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(keys);
460	>	if (inversionType) {
461	>	return inversionType;
462	>	} else {
463	>
464	>	AtomType* atype1;
465	>	AtomType* atype2;
466	>	AtomType* atype3;
467	>	AtomType* atype4;
468	>	std::vector<std::string> at1key;
469	>	at1key.push_back(at1);
470	>	atype1 = atomTypeCont_.find(at1key);
471	>
472	>	std::vector<std::string> at2key;
473	>	at2key.push_back(at2);
474	>	atype2 = atomTypeCont_.find(at2key);
475	>
476	>	std::vector<std::string> at3key;
477	>	at3key.push_back(at3);
478	>	atype3 = atomTypeCont_.find(at3key);
479	>
480	>	std::vector<std::string> at4key;
481	>	at4key.push_back(at4);
482	>	atype4 = atomTypeCont_.find(at4key);
483
484	+	// query atom types for their chains of responsibility
485	+	std::vector<AtomType*> at1Chain = atype1->allYourBase();
486	+	std::vector<AtomType*> at2Chain = atype2->allYourBase();
487	+	std::vector<AtomType*> at3Chain = atype3->allYourBase();
488	+	std::vector<AtomType*> at4Chain = atype4->allYourBase();
489	+
490	+	std::vector<AtomType*>::iterator i;
491	+	std::vector<AtomType*>::iterator j;
492	+	std::vector<AtomType*>::iterator k;
493	+	std::vector<AtomType*>::iterator l;
494	+
495	+	int ii = 0;
496	+	int jj = 0;
497	+	int kk = 0;
498	+	int ll = 0;
499	+	int Iscore;
500	+	int JKLscore;
501
502	<	TorsionType* bestType = torsionTypeCont_.find(theKeys);
503	<	if (bestType) {
504	<	return bestType;
502	>	std::vector<tuple3<int, int, std::vector<std::string> > > foundInversions;
503	>
504	>	for (j = at2Chain.begin(); j != at2Chain.end(); ++j) {
505	>	kk = 0;
506	>	for (k = at3Chain.begin(); k != at3Chain.end(); ++k) {
507	>	ii = 0;
508	>	for (i = at1Chain.begin(); i != at1Chain.end(); ++i) {
509	>	ll = 0;
510	>	for (l = at4Chain.begin(); l != at4Chain.end(); ++l) {
511	>
512	>	Iscore = ii;
513	>	JKLscore = jj + kk + ll;
514	>
515	>	std::vector<std::string> myKeys;
516	>	myKeys.push_back((*i)->getName());
517	>	myKeys.push_back((*j)->getName());
518	>	myKeys.push_back((*k)->getName());
519	>	myKeys.push_back((*l)->getName());
520	>
521	>	InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(myKeys);
522	>	if (inversionType) {
523	>	foundInversions.push_back( make_tuple3(Iscore, JKLscore, myKeys) );
524	>	}
525	>	ll++;
526	>	}
527	>	ii++;
528	>	}
529	>	kk++;
530	>	}
531	>	jj++;
532	>	}
533	>
534	>	if (!foundInversions.empty()) {
535	>	std::sort(foundInversions.begin(), foundInversions.end());
536	>	std::vector<std::string> theKeys = foundInversions[0].third;
537	>
538	>	InversionType* bestType = inversionTypeCont_.permutedFindSkippingFirstElement(theKeys);
539	>	return bestType;
540		} else {
541		//if no exact match found, try wild card match
542	<	return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
542	>	return inversionTypeCont_.find(keys, wildCardAtomTypeName_);
543		}
544		}
545		}
546	<
546	>
547		NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
548	+
549		std::vector<std::string> keys;
550		keys.push_back(at1);
551		keys.push_back(at2);
#	Line 356 | Line 555 | namespace oopse {
555		if (nbiType) {
556		return nbiType;
557		} else {
558	<	//if no exact match found, try wild card match
559	<	return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
560	<	}
558	>	AtomType* atype1;
559	>	AtomType* atype2;
560	>	std::vector<std::string> at1key;
561	>	at1key.push_back(at1);
562	>	atype1 = atomTypeCont_.find(at1key);
563	>
564	>	std::vector<std::string> at2key;
565	>	at2key.push_back(at2);
566	>	atype2 = atomTypeCont_.find(at2key);
567	>
568	>	// query atom types for their chains of responsibility
569	>	std::vector<AtomType*> at1Chain = atype1->allYourBase();
570	>	std::vector<AtomType*> at2Chain = atype2->allYourBase();
571	>
572	>	std::vector<AtomType*>::iterator i;
573	>	std::vector<AtomType*>::iterator j;
574	>
575	>	int ii = 0;
576	>	int jj = 0;
577	>	int nbiTypeScore;
578	>
579	>	std::vector<std::pair<int, std::vector<std::string> > > foundNBI;
580	>
581	>	for (i = at1Chain.begin(); i != at1Chain.end(); ++i) {
582	>	jj = 0;
583	>	for (j = at2Chain.begin(); j != at2Chain.end(); ++j) {
584	>
585	>	nbiTypeScore = ii + jj;
586	>
587	>	std::vector<std::string> myKeys;
588	>	myKeys.push_back((*i)->getName());
589	>	myKeys.push_back((*j)->getName());
590	>
591	>	NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(myKeys);
592	>	if (nbiType) {
593	>	foundNBI.push_back(std::make_pair(nbiTypeScore, myKeys));
594	>	}
595	>	jj++;
596	>	}
597	>	ii++;
598	>	}
599	>
600	>
601	>	if (!foundNBI.empty()) {
602	>	// sort the foundNBI by the score:
603	>	std::sort(foundNBI.begin(), foundNBI.end());
604	>	std::vector<std::string> theKeys = foundNBI[0].second;
605	>
606	>	NonBondedInteractionType* bestType = nonBondedInteractionTypeCont_.find(theKeys);
607	>	return bestType;
608	>	} else {
609	>	//if no exact match found, try wild card match
610	>	return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
611	>	}
612	>	}
613		}
614
615		BondType* ForceField::getExactBondType(const std::string &at1,
#	Line 390 | Line 641 | namespace oopse {
641		keys.push_back(at4);
642		return torsionTypeCont_.find(keys);
643		}
644	<
644	>
645	>	InversionType* ForceField::getExactInversionType(const std::string &at1,
646	>	const std::string &at2,
647	>	const std::string &at3,
648	>	const std::string &at4){
649	>	std::vector<std::string> keys;
650	>	keys.push_back(at1);
651	>	keys.push_back(at2);
652	>	keys.push_back(at3);
653	>	keys.push_back(at4);
654	>	return inversionTypeCont_.find(keys);
655	>	}
656	>
657		NonBondedInteractionType* ForceField::getExactNonBondedInteractionType(const std::string &at1, const std::string &at2){
658		std::vector<std::string> keys;
659		keys.push_back(at1);
660		keys.push_back(at2);
661		return nonBondedInteractionTypeCont_.find(keys);
662		}
663	+
664
401	–
665		bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
666		std::vector<std::string> keys;
667		keys.push_back(at);
668	+	atypeIdentToName[atomType->getIdent()] = at;
669		return atomTypeCont_.add(keys, atomType);
670		}
671
672	+	bool ForceField::replaceAtomType(const std::string &at, AtomType* atomType) {
673	+	std::vector<std::string> keys;
674	+	keys.push_back(at);
675	+	atypeIdentToName[atomType->getIdent()] = at;
676	+	return atomTypeCont_.replace(keys, atomType);
677	+	}
678	+
679		bool ForceField::addBondType(const std::string &at1, const std::string &at2,
680		BondType* bondType) {
681		std::vector<std::string> keys;
#	Line 435 | Line 706 | namespace oopse {
706		return torsionTypeCont_.add(keys, torsionType);
707		}
708
709	+	bool ForceField::addInversionType(const std::string &at1,
710	+	const std::string &at2,
711	+	const std::string &at3,
712	+	const std::string &at4,
713	+	InversionType* inversionType) {
714	+	std::vector<std::string> keys;
715	+	keys.push_back(at1);
716	+	keys.push_back(at2);
717	+	keys.push_back(at3);
718	+	keys.push_back(at4);
719	+	return inversionTypeCont_.add(keys, inversionType);
720	+	}
721	+
722		bool ForceField::addNonBondedInteractionType(const std::string &at1,
723		const std::string &at2,
724		NonBondedInteractionType* nbiType) {
#	Line 445 | Line 729 | namespace oopse {
729		}
730
731		RealType ForceField::getRcutFromAtomType(AtomType* at) {
732	<	/**@todo */
449	<	GenericData* data;
450	<	RealType rcut = 0.0;
732	>	RealType rcut(0.0);
733
734	<	if (at->isLennardJones()) {
735	<	data = at->getPropertyByName("LennardJones");
736	<	if (data != NULL) {
455	<	LJParamGenericData* ljData = dynamic_cast<LJParamGenericData*>(data);
456	<
457	<	if (ljData != NULL) {
458	<	LJParam ljParam = ljData->getData();
459	<
460	<	//by default use 2.5*sigma as cutoff radius
461	<	rcut = 2.5 * ljParam.sigma;
462	<
463	<	} else {
464	<	sprintf( painCave.errMsg,
465	<	"Can not cast GenericData to LJParam\n");
466	<	painCave.severity = OOPSE_ERROR;
467	<	painCave.isFatal = 1;
468	<	simError();
469	<	}
470	<	} else {
471	<	sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
472	<	painCave.severity = OOPSE_ERROR;
473	<	painCave.isFatal = 1;
474	<	simError();
475	<	}
734	>	LennardJonesAdapter lja = LennardJonesAdapter(at);
735	>	if (lja.isLennardJones()) {
736	>	rcut = 2.5 * lja.getSigma();
737		}
738	+	EAMAdapter ea = EAMAdapter(at);
739	+	if (ea.isEAM()) {
740	+	rcut = max(rcut, ea.getRcut());
741	+	}
742	+	SuttonChenAdapter sca = SuttonChenAdapter(at);
743	+	if (sca.isSuttonChen()) {
744	+	rcut = max(rcut, 2.0 * sca.getAlpha());
745	+	}
746	+	GayBerneAdapter gba = GayBerneAdapter(at);
747	+	if (gba.isGayBerne()) {
748	+	rcut = max(rcut, 2.5 * sqrt(2.0) * max(gba.getD(), gba.getL()));
749	+	}
750	+	StickyAdapter sa = StickyAdapter(at);
751	+	if (sa.isSticky()) {
752	+	rcut = max(rcut, max(sa.getRu(), sa.getRup()));
753	+	}
754	+
755		return rcut;
756		}
757
#	Line 499 | Line 777 | namespace oopse {
777		"\tHave you tried setting the FORCE_PARAM_PATH environment "
778		"variable?\n",
779		forceFieldFilename.c_str() );
780	<	painCave.severity = OOPSE_ERROR;
780	>	painCave.severity = OPENMD_ERROR;
781		painCave.isFatal = 1;
782		simError();
783		}
#	Line 507 | Line 785 | namespace oopse {
785		return ffStream;
786		}
787
788	<	void ForceField::setFortranForceOptions(){
511	<	ForceOptions theseFortranOptions;
512	<	forceFieldOptions_.makeFortranOptions(theseFortranOptions);
513	<	setfForceOptions(&theseFortranOptions);
514	<	}
515	<	} //end namespace oopse
788	>	} //end namespace OpenMD

Comparing:
trunk/src/UseTheForce/ForceField.cpp (property svn:keywords), Revision 1269 by gezelter, Tue Jul 1 13:28:23 2008 UTC vs.
trunk/src/brains/ForceField.cpp (property svn:keywords), Revision 1880 by gezelter, Mon Jun 17 18:28:30 2013 UTC

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

Diff Legend

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