51 |
|
|
52 |
|
P2OrderParameter::P2OrderParameter(SimInfo* info, const string& filename, |
53 |
|
const string& sele1) |
54 |
< |
: StaticAnalyser(info, filename), doVect_(true), |
55 |
< |
selectionScript1_(sele1), evaluator1_(info), |
56 |
< |
evaluator2_(info), seleMan1_(info), seleMan2_(info) { |
54 |
> |
: StaticAnalyser(info, filename), doVect_(true), doOffset_(false), |
55 |
> |
selectionScript1_(sele1), evaluator1_(info), |
56 |
> |
evaluator2_(info), seleMan1_(info), seleMan2_(info) { |
57 |
|
|
58 |
|
setOutputName(getPrefix(filename) + ".p2"); |
59 |
|
|
62 |
|
|
63 |
|
P2OrderParameter::P2OrderParameter(SimInfo* info, const string& filename, |
64 |
|
const string& sele1, const string& sele2) |
65 |
< |
: StaticAnalyser(info, filename), doVect_(false), |
66 |
< |
selectionScript1_(sele1), selectionScript2_(sele2), evaluator1_(info), |
67 |
< |
evaluator2_(info), seleMan1_(info), seleMan2_(info) { |
65 |
> |
: StaticAnalyser(info, filename), doVect_(false), doOffset_(false), |
66 |
> |
selectionScript1_(sele1), selectionScript2_(sele2), evaluator1_(info), |
67 |
> |
evaluator2_(info), seleMan1_(info), seleMan2_(info) { |
68 |
|
|
69 |
|
setOutputName(getPrefix(filename) + ".p2"); |
70 |
|
|
72 |
|
evaluator2_.loadScriptString(sele2); |
73 |
|
} |
74 |
|
|
75 |
+ |
P2OrderParameter::P2OrderParameter(SimInfo* info, const string& filename, |
76 |
+ |
const string& sele1, int seleOffset) |
77 |
+ |
: StaticAnalyser(info, filename), doVect_(false), doOffset_(true), |
78 |
+ |
seleOffset_(seleOffset), selectionScript1_(sele1), evaluator1_(info), |
79 |
+ |
evaluator2_(info), seleMan1_(info), seleMan2_(info) { |
80 |
+ |
|
81 |
+ |
setOutputName(getPrefix(filename) + ".p2"); |
82 |
+ |
|
83 |
+ |
evaluator1_.loadScriptString(sele1); |
84 |
+ |
} |
85 |
+ |
|
86 |
|
void P2OrderParameter::process() { |
87 |
|
Molecule* mol; |
88 |
|
RigidBody* rb; |
113 |
|
Mat3x3d orderTensor(0.0); |
114 |
|
vecCount = 0; |
115 |
|
|
116 |
< |
if (evaluator1_.isDynamic()) |
106 |
< |
seleMan1_.setSelectionSet(evaluator1_.evaluate()); |
116 |
> |
seleMan1_.setSelectionSet(evaluator1_.evaluate()); |
117 |
|
|
118 |
|
if (doVect_) { |
119 |
|
|
131 |
|
|
132 |
|
} else { |
133 |
|
|
134 |
< |
if (evaluator2_.isDynamic()) |
125 |
< |
seleMan2_.setSelectionSet(evaluator2_.evaluate()); |
126 |
< |
|
127 |
< |
if (seleMan1_.getSelectionCount() != seleMan2_.getSelectionCount() ) { |
128 |
< |
sprintf( painCave.errMsg, |
129 |
< |
"In frame %d, the number of selected StuntDoubles are not\n" |
130 |
< |
"\tthe same in --sele1 and sele2\n", i); |
131 |
< |
painCave.severity = OPENMD_INFO; |
132 |
< |
painCave.isFatal = 0; |
133 |
< |
simError(); |
134 |
< |
} |
135 |
< |
|
136 |
< |
for (sd1 = seleMan1_.beginSelected(ii), |
137 |
< |
sd2 = seleMan2_.beginSelected(jj); |
138 |
< |
sd1 != NULL && sd2 != NULL; |
139 |
< |
sd1 = seleMan1_.nextSelected(ii), |
140 |
< |
sd2 = seleMan2_.nextSelected(jj)) { |
134 |
> |
if (doOffset_) { |
135 |
|
|
136 |
< |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
136 |
> |
for (sd1 = seleMan1_.beginSelected(ii); sd1 != NULL; |
137 |
> |
sd1 = seleMan1_.nextSelected(ii)) { |
138 |
|
|
139 |
< |
if (usePeriodicBoundaryConditions_) |
140 |
< |
currentSnapshot_->wrapVector(vec); |
139 |
> |
// This will require careful rewriting if StaticProps is |
140 |
> |
// ever parallelized. For an example, see |
141 |
> |
// Thermo::getTaggedAtomPairDistance |
142 |
|
|
143 |
< |
vec.normalize(); |
143 |
> |
int sd2Index = sd1->getGlobalIndex() + seleOffset_; |
144 |
> |
sd2 = info_->getIOIndexToIntegrableObject(sd2Index); |
145 |
> |
|
146 |
> |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
147 |
> |
|
148 |
> |
if (usePeriodicBoundaryConditions_) |
149 |
> |
currentSnapshot_->wrapVector(vec); |
150 |
> |
|
151 |
> |
vec.normalize(); |
152 |
> |
|
153 |
> |
orderTensor +=outProduct(vec, vec); |
154 |
> |
vecCount++; |
155 |
> |
} |
156 |
> |
|
157 |
> |
orderTensor /= vecCount; |
158 |
> |
} else { |
159 |
> |
|
160 |
> |
seleMan2_.setSelectionSet(evaluator2_.evaluate()); |
161 |
> |
|
162 |
> |
if (seleMan1_.getSelectionCount() != seleMan2_.getSelectionCount() ) { |
163 |
> |
sprintf( painCave.errMsg, |
164 |
> |
"In frame %d, the number of selected StuntDoubles are\n" |
165 |
> |
"\tnot the same in --sele1 and sele2\n", i); |
166 |
> |
painCave.severity = OPENMD_INFO; |
167 |
> |
painCave.isFatal = 0; |
168 |
> |
simError(); |
169 |
> |
} |
170 |
> |
|
171 |
> |
for (sd1 = seleMan1_.beginSelected(ii), |
172 |
> |
sd2 = seleMan2_.beginSelected(jj); |
173 |
> |
sd1 != NULL && sd2 != NULL; |
174 |
> |
sd1 = seleMan1_.nextSelected(ii), |
175 |
> |
sd2 = seleMan2_.nextSelected(jj)) { |
176 |
> |
|
177 |
> |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
178 |
> |
|
179 |
> |
if (usePeriodicBoundaryConditions_) |
180 |
> |
currentSnapshot_->wrapVector(vec); |
181 |
|
|
182 |
< |
orderTensor +=outProduct(vec, vec); |
183 |
< |
vecCount++; |
182 |
> |
vec.normalize(); |
183 |
> |
|
184 |
> |
orderTensor +=outProduct(vec, vec); |
185 |
> |
vecCount++; |
186 |
> |
} |
187 |
> |
|
188 |
> |
orderTensor /= vecCount; |
189 |
|
} |
152 |
– |
|
153 |
– |
orderTensor /= vecCount; |
190 |
|
} |
191 |
|
|
192 |
|
if (vecCount == 0) { |
237 |
|
angle = angle/(vecCount*NumericConstant::PI)*180.0; |
238 |
|
|
239 |
|
} else { |
240 |
< |
for (sd1 = seleMan1_.beginSelected(ii), |
241 |
< |
sd2 = seleMan2_.beginSelected(jj); |
242 |
< |
sd1 != NULL && sd2 != NULL; |
243 |
< |
sd1 = seleMan1_.nextSelected(ii), |
244 |
< |
sd2 = seleMan2_.nextSelected(jj)) { |
245 |
< |
|
246 |
< |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
247 |
< |
if (usePeriodicBoundaryConditions_) |
248 |
< |
currentSnapshot_->wrapVector(vec); |
249 |
< |
vec.normalize(); |
250 |
< |
angle += acos(dot(vec, director)) ; |
251 |
< |
vecCount++; |
240 |
> |
if (doOffset_) { |
241 |
> |
|
242 |
> |
for (sd1 = seleMan1_.beginSelected(ii); sd1 != NULL; |
243 |
> |
sd1 = seleMan1_.nextSelected(ii)) { |
244 |
> |
|
245 |
> |
// This will require careful rewriting if StaticProps is |
246 |
> |
// ever parallelized. For an example, see |
247 |
> |
// Thermo::getTaggedAtomPairDistance |
248 |
> |
|
249 |
> |
int sd2Index = sd1->getGlobalIndex() + seleOffset_; |
250 |
> |
sd2 = info_->getIOIndexToIntegrableObject(sd2Index); |
251 |
> |
|
252 |
> |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
253 |
> |
if (usePeriodicBoundaryConditions_) |
254 |
> |
currentSnapshot_->wrapVector(vec); |
255 |
> |
vec.normalize(); |
256 |
> |
angle += acos(dot(vec, director)) ; |
257 |
> |
vecCount++; |
258 |
> |
} |
259 |
> |
angle = angle / (vecCount * NumericConstant::PI) * 180.0; |
260 |
> |
|
261 |
> |
} else { |
262 |
> |
|
263 |
> |
for (sd1 = seleMan1_.beginSelected(ii), |
264 |
> |
sd2 = seleMan2_.beginSelected(jj); |
265 |
> |
sd1 != NULL && sd2 != NULL; |
266 |
> |
sd1 = seleMan1_.nextSelected(ii), |
267 |
> |
sd2 = seleMan2_.nextSelected(jj)) { |
268 |
> |
|
269 |
> |
Vector3d vec = sd1->getPos() - sd2->getPos(); |
270 |
> |
if (usePeriodicBoundaryConditions_) |
271 |
> |
currentSnapshot_->wrapVector(vec); |
272 |
> |
vec.normalize(); |
273 |
> |
angle += acos(dot(vec, director)) ; |
274 |
> |
vecCount++; |
275 |
> |
} |
276 |
> |
angle = angle / (vecCount * NumericConstant::PI) * 180.0; |
277 |
|
} |
217 |
– |
angle = angle / (vecCount * NumericConstant::PI) * 180.0; |
278 |
|
} |
279 |
|
|
280 |
|
OrderParam param; |