ViewVC Help

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

root/OpenMD/trunk/src/applications/staticProps/GofRAngle.cpp

Comparing trunk/src/applications/staticProps/GofRAngle.cpp (file contents):
Revision 963 by tim, Wed May 17 21:51:42 2006 UTC vs.
Revision 1968 by gezelter, Thu Feb 20 16:27:30 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 <algorithm>
44		#include <fstream>
45		#include "applications/staticProps/GofRAngle.hpp"
46	+	#include "primitives/Atom.hpp"
47	+	#include "types/MultipoleAdapter.hpp"
48		#include "utils/simError.h"
49
50	<	namespace oopse {
50	>	namespace OpenMD {
51
52		GofRAngle::GofRAngle(SimInfo* info, const std::string& filename, const std::string& sele1,
53		const std::string& sele2, RealType len, int nrbins, int nangleBins)
54		: RadialDistrFunc(info, filename, sele1, sele2), len_(len), nRBins_(nrbins), nAngleBins_(nangleBins){
55
56	<	deltaR_ = len_ /nRBins_;
57	<	deltaCosAngle_ = 2.0 / nAngleBins_;
55	<
56	>	deltaR_ = len_ /(double) nRBins_;
57	>	deltaCosAngle_ = 2.0 / (double)nAngleBins_;
58		histogram_.resize(nRBins_);
59		avgGofr_.resize(nRBins_);
60		for (int i = 0 ; i < nRBins_; ++i) {
#	Line 63 | Line 65 | namespace oopse {
65
66
67		void GofRAngle::preProcess() {
68	<
67	<	for (int i = 0; i < avgGofr_.size(); ++i) {
68	>	for (unsigned int i = 0; i < avgGofr_.size(); ++i) {
69		std::fill(avgGofr_[i].begin(), avgGofr_[i].end(), 0);
70		}
71		}
72
73	<	void GofRAngle::initalizeHistogram() {
73	>	void GofRAngle::initializeHistogram() {
74		npairs_ = 0;
75	<	for (int i = 0; i < histogram_.size(); ++i)
75	>	for (unsigned int i = 0; i < histogram_.size(); ++i){
76		std::fill(histogram_[i].begin(), histogram_[i].end(), 0);
77	+	}
78		}
79
78	–
80		void GofRAngle::processHistogram() {
80	–
81		int nPairs = getNPairs();
82		RealType volume = info_->getSnapshotManager()->getCurrentSnapshot()->getVolume();
83		RealType pairDensity = nPairs /volume;
84		RealType pairConstant = ( 4.0 * NumericConstant::PI * pairDensity ) / 3.0;
85
86	<	for(int i = 0 ; i < histogram_.size(); ++i){
86	>	for(unsigned int i = 0 ; i < histogram_.size(); ++i){
87
88		RealType rLower = i * deltaR_;
89		RealType rUpper = rLower + deltaR_;
90		RealType volSlice = ( rUpper * rUpper * rUpper ) - ( rLower * rLower * rLower );
91		RealType nIdeal = volSlice * pairConstant;
92
93	<	for (int j = 0; j < histogram_[i].size(); ++j){
93	>	for (unsigned int j = 0; j < histogram_[i].size(); ++j){
94		avgGofr_[i][j] += histogram_[i][j] / nIdeal;
95		}
96		}
#	Line 102 | Line 102 | namespace oopse {
102		if (sd1 == sd2) {
103		return;
104		}
105	–
105		Vector3d pos1 = sd1->getPos();
106		Vector3d pos2 = sd2->getPos();
107		Vector3d r12 = pos2 - pos1;
108	<	currentSnapshot_->wrapVector(r12);
108	>	if (usePeriodicBoundaryConditions_)
109	>	currentSnapshot_->wrapVector(r12);
110
111		RealType distance = r12.length();
112	<	int whichRBin = distance / deltaR_;
112	>	int whichRBin = int(distance / deltaR_);
113
114		if (distance <= len_) {
115	+
116		RealType cosAngle = evaluateAngle(sd1, sd2);
117		RealType halfBin = (nAngleBins_ - 1) * 0.5;
118	<	int whichThetaBin = halfBin * (cosAngle + 1.0);
118	>	int whichThetaBin = int(halfBin * (cosAngle + 1.0));
119		++histogram_[whichRBin][whichThetaBin];
120
121		++npairs_;
#	Line 129 | Line 130 | namespace oopse {
130		rdfStream << "selection2: (" << selectionScript2_ << ")\n";
131		rdfStream << "#nRBins = " << nRBins_ << "\t maxLen = " << len_ << "deltaR = " << deltaR_ <<"\n";
132		rdfStream << "#nAngleBins =" << nAngleBins_ << "deltaCosAngle = " << deltaCosAngle_ << "\n";
133	<	for (int i = 0; i < avgGofr_.size(); ++i) {
134	<	RealType r = deltaR_ * (i + 0.5);
133	>	for (unsigned int i = 0; i < avgGofr_.size(); ++i) {
134	>	// RealType r = deltaR_ * (i + 0.5);
135
136	<	for(int j = 0; j < avgGofr_[i].size(); ++j) {
137	<	RealType cosAngle = -1.0 + (j + 0.5)*deltaCosAngle_;
136	>	for(unsigned int j = 0; j < avgGofr_[i].size(); ++j) {
137	>	// RealType cosAngle = -1.0 + (j + 0.5)*deltaCosAngle_;
138		rdfStream << avgGofr_[i][j]/nProcessed_ << "\t";
139		}
140
#	Line 153 | Line 154 | namespace oopse {
154		Vector3d pos1 = sd1->getPos();
155		Vector3d pos2 = sd2->getPos();
156		Vector3d r12 = pos2 - pos1;
157	<	currentSnapshot_->wrapVector(r12);
157	>
158	>	if (usePeriodicBoundaryConditions_)
159	>	currentSnapshot_->wrapVector(r12);
160	>
161		r12.normalize();
162	<	Vector3d dipole = sd1->getElectroFrame().getColumn(2);
163	<	dipole.normalize();
164	<	return dot(r12, dipole);
162	>
163	>	Vector3d vec;
164	>
165	>	if (sd1->isAtom()) {
166	>	AtomType* atype1 = static_cast<Atom*>(sd1)->getAtomType();
167	>	MultipoleAdapter ma1 = MultipoleAdapter(atype1);
168	>
169	>	if (ma1.isDipole() )
170	>	vec = sd1->getDipole();
171	>	else
172	>	vec = sd1->getA().transpose() * V3Z;
173	>	} else {
174	>	vec = sd1->getA().transpose() * V3Z;
175	>	}
176	>
177	>	vec.normalize();
178	>
179	>	return dot(r12, vec);
180		}
181
182		RealType GofROmega::evaluateAngle(StuntDouble* sd1, StuntDouble* sd2) {
183	<	Vector3d v1 = sd1->getElectroFrame().getColumn(2);
184	<	Vector3d v2 = sd2->getElectroFrame().getColumn(2);
183	>	Vector3d v1, v2;
184	>
185	>	if (sd1->isAtom()){
186	>	AtomType* atype1 = static_cast<Atom*>(sd1)->getAtomType();
187	>	MultipoleAdapter ma1 = MultipoleAdapter(atype1);
188	>	if (ma1.isDipole() )
189	>	v1 = sd1->getDipole();
190	>	else
191	>	v1 = sd1->getA().transpose() * V3Z;
192	>	} else {
193	>	v1 = sd1->getA().transpose() * V3Z;
194	>	}
195	>
196	>	if (sd2->isAtom()) {
197	>	AtomType* atype2 = static_cast<Atom*>(sd2)->getAtomType();
198	>	MultipoleAdapter ma2 = MultipoleAdapter(atype2);
199	>
200	>	if (ma2.isDipole() )
201	>	v2 = sd2->getDipole();
202	>	else
203	>	v2 = sd2->getA().transpose() * V3Z;
204	>	} else {
205	>	v2 = sd2->getA().transpose() * V3Z;
206	>	}
207	>
208		v1.normalize();
209		v2.normalize();
210		return dot(v1, v2);

Comparing trunk/src/applications/staticProps/GofRAngle.cpp (property svn:keywords):
Revision 963 by tim, Wed May 17 21:51:42 2006 UTC vs.
Revision 1968 by gezelter, Thu Feb 20 16:27:30 2014 UTC

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

Diff Legend

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