source: src/analysis_correlation.cpp@ f2a1d3

Action_Thermostats Add_AtomRandomPerturbation Add_FitFragmentPartialChargesAction Add_RotateAroundBondAction Add_SelectAtomByNameAction Added_ParseSaveFragmentResults AddingActions_SaveParseParticleParameters Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_ParticleName_to_Atom Adding_StructOpt_integration_tests AtomFragments Automaking_mpqc_open AutomationFragmentation_failures Candidate_v1.5.4 Candidate_v1.6.0 Candidate_v1.6.1 Candidate_v1.7.0 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator CombiningParticlePotentialParsing Combining_Subpackages Debian_Package_split Debian_package_split_molecuildergui_only Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph EmpiricalPotential_contain_HomologyGraph_documentation Enable_parallel_make_install Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_BoundInBox_CenterInBox_MoleculeActions Fix_ChargeSampling_PBC Fix_ChronosMutex Fix_FitPartialCharges Fix_FitPotential_needs_atomicnumbers Fix_ForceAnnealing Fix_IndependentFragmentGrids Fix_ParseParticles Fix_ParseParticles_split_forward_backward_Actions Fix_PopActions Fix_QtFragmentList_sorted_selection Fix_Restrictedkeyset_FragmentMolecule Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns Fix_fitting_potentials Fixes ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion FragmentAction_writes_AtomFragments FragmentMolecule_checks_bonddegrees GeometryObjects Gui_Fixes Gui_displays_atomic_force_velocity ImplicitCharges IndependentFragmentGrids IndependentFragmentGrids_IndividualZeroInstances IndependentFragmentGrids_IntegrationTest IndependentFragmentGrids_Sole_NN_Calculation JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool JobMarket_unresolvable_hostname_fix MoreRobust_FragmentAutomation ODR_violation_mpqc_open PartialCharges_OrthogonalSummation PdbParser_setsAtomName PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks Rewrite_FitPartialCharges RotateToPrincipalAxisSystem_UndoRedo SaturateAtoms_findBestMatching SaturateAtoms_singleDegree StoppableMakroAction Subpackage_CodePatterns Subpackage_JobMarket Subpackage_LinearAlgebra Subpackage_levmar Subpackage_mpqc_open Subpackage_vmg Switchable_LogView ThirdParty_MPQC_rebuilt_buildsystem TrajectoryDependenant_MaxOrder TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps TremoloParser_setsAtomName Ubuntu_1604_changes stable
Last change on this file since f2a1d3 was 7fd416, checked in by Frederik Heber <heber@…>, 15 years ago

FIX: CorrelationToSurface() was broken.

  • DOCU: CorrelationToSurface() is more verbose on empty molecules, molecule::doCountAtoms() is verbose on naming atoms only for high verbosity levels
  • MEMFIX: ParseCommandLineOptions() - case 'CP' did not set counter to zero for re-setting active flag.
  • BUGFIX: molecule::DeterminePeriodicCenter() - InverseMatrix() was called with cell_size instead of full matrix.
  • rewritten MoleculeListClass::DissectMoleculeIntoConnectedSubgraphs() a bit:
    • molname and contained atoms are given
    • no more stupid map, atoms are directly transfered from Leaf to molecules[]
    • mol that contained all was not destroyed after use
  • Property mode set to 100644
File size: 20.4 KB
Line 
1/*
2 * analysis.cpp
3 *
4 * Created on: Oct 13, 2009
5 * Author: heber
6 */
7
8#include <iostream>
9
10#include "analysis_correlation.hpp"
11#include "element.hpp"
12#include "info.hpp"
13#include "log.hpp"
14#include "molecule.hpp"
15#include "tesselation.hpp"
16#include "tesselationhelpers.hpp"
17#include "triangleintersectionlist.hpp"
18#include "vector.hpp"
19#include "verbose.hpp"
20#include "World.hpp"
21
22
23/** Calculates the pair correlation between given elements.
24 * Note given element order is unimportant (i.e. g(Si, O) === g(O, Si))
25 * \param *out output stream for debugging
26 * \param *molecules list of molecules structure
27 * \param *type1 first element or NULL (if any element)
28 * \param *type2 second element or NULL (if any element)
29 * \return Map of doubles with values the pair of the two atoms.
30 */
31PairCorrelationMap *PairCorrelation(MoleculeListClass * const &molecules, const element * const type1, const element * const type2 )
32{
33 Info FunctionInfo(__func__);
34 PairCorrelationMap *outmap = NULL;
35 double distance = 0.;
36
37 if (molecules->ListOfMolecules.empty()) {
38 DoeLog(1) && (eLog()<< Verbose(1) <<"No molecule given." << endl);
39 return outmap;
40 }
41 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
42 (*MolWalker)->doCountAtoms();
43 outmap = new PairCorrelationMap;
44 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++){
45 if ((*MolWalker)->ActiveFlag) {
46 DoeLog(2) && (eLog()<< Verbose(2) << "Current molecule is " << *MolWalker << "." << endl);
47 eLog() << Verbose(2) << "Current molecule is " << *MolWalker << "." << endl;
48 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
49 DoLog(3) && (Log() << Verbose(3) << "Current atom is " << **iter << "." << endl);
50 if ((type1 == NULL) || ((*iter)->type == type1)) {
51 for (MoleculeList::const_iterator MolOtherWalker = MolWalker; MolOtherWalker != molecules->ListOfMolecules.end(); MolOtherWalker++){
52 if ((*MolOtherWalker)->ActiveFlag) {
53 DoLog(2) && (Log() << Verbose(2) << "Current other molecule is " << *MolOtherWalker << "." << endl);
54 for (molecule::const_iterator runner = (*MolOtherWalker)->begin(); runner != (*MolOtherWalker)->end(); ++runner) {
55 DoLog(3) && (Log() << Verbose(3) << "Current otheratom is " << **runner << "." << endl);
56 if ((*iter)->getId() < (*runner)->getId()){
57 if ((type2 == NULL) || ((*runner)->type == type2)) {
58 distance = (*iter)->node->PeriodicDistance(*(*runner)->node, World::getInstance().getDomain());
59 //Log() << Verbose(1) <<"Inserting " << *(*iter) << " and " << *(*runner) << endl;
60 outmap->insert ( pair<double, pair <atom *, atom*> > (distance, pair<atom *, atom*> ((*iter), (*runner)) ) );
61 }
62 }
63 }
64 }
65 }
66 }
67 }
68 }
69 }
70 return outmap;
71};
72
73/** Calculates the pair correlation between given elements.
74 * Note given element order is unimportant (i.e. g(Si, O) === g(O, Si))
75 * \param *out output stream for debugging
76 * \param *molecules list of molecules structure
77 * \param *type1 first element or NULL (if any element)
78 * \param *type2 second element or NULL (if any element)
79 * \param ranges[NDIM] interval boundaries for the periodic images to scan also
80 * \return Map of doubles with values the pair of the two atoms.
81 */
82PairCorrelationMap *PeriodicPairCorrelation(MoleculeListClass * const &molecules, const element * const type1, const element * const type2, const int ranges[NDIM] )
83{
84 Info FunctionInfo(__func__);
85 PairCorrelationMap *outmap = NULL;
86 double distance = 0.;
87 int n[NDIM];
88 Vector checkX;
89 Vector periodicX;
90 int Othern[NDIM];
91 Vector checkOtherX;
92 Vector periodicOtherX;
93
94 if (molecules->ListOfMolecules.empty()) {
95 DoeLog(1) && (eLog()<< Verbose(1) <<"No molecule given." << endl);
96 return outmap;
97 }
98 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
99 (*MolWalker)->doCountAtoms();
100 outmap = new PairCorrelationMap;
101 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
102 if ((*MolWalker)->ActiveFlag) {
103 double * FullMatrix = ReturnFullMatrixforSymmetric(World::getInstance().getDomain());
104 double * FullInverseMatrix = InverseMatrix(FullMatrix);
105 DoeLog(2) && (eLog()<< Verbose(2) << "Current molecule is " << *MolWalker << "." << endl);
106 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
107 DoLog(3) && (Log() << Verbose(3) << "Current atom is " << **iter << "." << endl);
108 if ((type1 == NULL) || ((*iter)->type == type1)) {
109 periodicX = *(*iter)->node;
110 periodicX.MatrixMultiplication(FullInverseMatrix); // x now in [0,1)^3
111 // go through every range in xyz and get distance
112 for (n[0]=-ranges[0]; n[0] <= ranges[0]; n[0]++)
113 for (n[1]=-ranges[1]; n[1] <= ranges[1]; n[1]++)
114 for (n[2]=-ranges[2]; n[2] <= ranges[2]; n[2]++) {
115 checkX = Vector(n[0], n[1], n[2]) + periodicX;
116 checkX.MatrixMultiplication(FullMatrix);
117 for (MoleculeList::const_iterator MolOtherWalker = MolWalker; MolOtherWalker != molecules->ListOfMolecules.end(); MolOtherWalker++)
118 if ((*MolOtherWalker)->ActiveFlag) {
119 DoLog(2) && (Log() << Verbose(2) << "Current other molecule is " << *MolOtherWalker << "." << endl);
120 for (molecule::const_iterator runner = (*MolOtherWalker)->begin(); runner != (*MolOtherWalker)->end(); ++runner) {
121 DoLog(3) && (Log() << Verbose(3) << "Current otheratom is " << **runner << "." << endl);
122 if ((*iter)->nr < (*runner)->nr)
123 if ((type2 == NULL) || ((*runner)->type == type2)) {
124 periodicOtherX = *(*runner)->node;
125 periodicOtherX.MatrixMultiplication(FullInverseMatrix); // x now in [0,1)^3
126 // go through every range in xyz and get distance
127 for (Othern[0]=-ranges[0]; Othern[0] <= ranges[0]; Othern[0]++)
128 for (Othern[1]=-ranges[1]; Othern[1] <= ranges[1]; Othern[1]++)
129 for (Othern[2]=-ranges[2]; Othern[2] <= ranges[2]; Othern[2]++) {
130 checkOtherX = Vector(Othern[0], Othern[1], Othern[2]) + periodicOtherX;
131 checkOtherX.MatrixMultiplication(FullMatrix);
132 distance = checkX.distance(checkOtherX);
133 //Log() << Verbose(1) <<"Inserting " << *(*iter) << " and " << *(*runner) << endl;
134 outmap->insert ( pair<double, pair <atom *, atom*> > (distance, pair<atom *, atom*> ((*iter), (*runner)) ) );
135 }
136 }
137 }
138 }
139 }
140 }
141 }
142 delete[](FullMatrix);
143 delete[](FullInverseMatrix);
144 }
145
146 return outmap;
147};
148
149/** Calculates the distance (pair) correlation between a given element and a point.
150 * \param *out output stream for debugging
151 * \param *molecules list of molecules structure
152 * \param *type element or NULL (if any element)
153 * \param *point vector to the correlation point
154 * \return Map of dobules with values as pairs of atom and the vector
155 */
156CorrelationToPointMap *CorrelationToPoint(MoleculeListClass * const &molecules, const element * const type, const Vector *point )
157{
158 Info FunctionInfo(__func__);
159 CorrelationToPointMap *outmap = NULL;
160 double distance = 0.;
161
162 if (molecules->ListOfMolecules.empty()) {
163 DoLog(1) && (Log() << Verbose(1) <<"No molecule given." << endl);
164 return outmap;
165 }
166 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
167 (*MolWalker)->doCountAtoms();
168 outmap = new CorrelationToPointMap;
169 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
170 if ((*MolWalker)->ActiveFlag) {
171 DoLog(2) && (Log() << Verbose(2) << "Current molecule is " << *MolWalker << "." << endl);
172 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
173 DoLog(3) && (Log() << Verbose(3) << "Current atom is " << **iter << "." << endl);
174 if ((type == NULL) || ((*iter)->type == type)) {
175 distance = (*iter)->node->PeriodicDistance(*point, World::getInstance().getDomain());
176 DoLog(4) && (Log() << Verbose(4) << "Current distance is " << distance << "." << endl);
177 outmap->insert ( pair<double, pair<atom *, const Vector*> >(distance, pair<atom *, const Vector*> ((*iter), point) ) );
178 }
179 }
180 }
181
182 return outmap;
183};
184
185/** Calculates the distance (pair) correlation between a given element, all its periodic images and a point.
186 * \param *out output stream for debugging
187 * \param *molecules list of molecules structure
188 * \param *type element or NULL (if any element)
189 * \param *point vector to the correlation point
190 * \param ranges[NDIM] interval boundaries for the periodic images to scan also
191 * \return Map of dobules with values as pairs of atom and the vector
192 */
193CorrelationToPointMap *PeriodicCorrelationToPoint(MoleculeListClass * const &molecules, const element * const type, const Vector *point, const int ranges[NDIM] )
194{
195 Info FunctionInfo(__func__);
196 CorrelationToPointMap *outmap = NULL;
197 double distance = 0.;
198 int n[NDIM];
199 Vector periodicX;
200 Vector checkX;
201
202 if (molecules->ListOfMolecules.empty()) {
203 DoLog(1) && (Log() << Verbose(1) <<"No molecule given." << endl);
204 return outmap;
205 }
206 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
207 (*MolWalker)->doCountAtoms();
208 outmap = new CorrelationToPointMap;
209 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
210 if ((*MolWalker)->ActiveFlag) {
211 double * FullMatrix = ReturnFullMatrixforSymmetric(World::getInstance().getDomain());
212 double * FullInverseMatrix = InverseMatrix(FullMatrix);
213 DoLog(2) && (Log() << Verbose(2) << "Current molecule is " << *MolWalker << "." << endl);
214 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
215 DoLog(3) && (Log() << Verbose(3) << "Current atom is " << **iter << "." << endl);
216 if ((type == NULL) || ((*iter)->type == type)) {
217 periodicX = *(*iter)->node;
218 periodicX.MatrixMultiplication(FullInverseMatrix); // x now in [0,1)^3
219 // go through every range in xyz and get distance
220 for (n[0]=-ranges[0]; n[0] <= ranges[0]; n[0]++)
221 for (n[1]=-ranges[1]; n[1] <= ranges[1]; n[1]++)
222 for (n[2]=-ranges[2]; n[2] <= ranges[2]; n[2]++) {
223 checkX = Vector(n[0], n[1], n[2]) + periodicX;
224 checkX.MatrixMultiplication(FullMatrix);
225 distance = checkX.distance(*point);
226 DoLog(4) && (Log() << Verbose(4) << "Current distance is " << distance << "." << endl);
227 outmap->insert ( pair<double, pair<atom *, const Vector*> >(distance, pair<atom *, const Vector*> (*iter, point) ) );
228 }
229 }
230 }
231 delete[](FullMatrix);
232 delete[](FullInverseMatrix);
233 }
234
235 return outmap;
236};
237
238/** Calculates the distance (pair) correlation between a given element and a surface.
239 * \param *out output stream for debugging
240 * \param *molecules list of molecules structure
241 * \param *type element or NULL (if any element)
242 * \param *Surface pointer to Tesselation class surface
243 * \param *LC LinkedCell structure to quickly find neighbouring atoms
244 * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest
245 */
246CorrelationToSurfaceMap *CorrelationToSurface(MoleculeListClass * const &molecules, const element * const type, const Tesselation * const Surface, const LinkedCell *LC )
247{
248 Info FunctionInfo(__func__);
249 CorrelationToSurfaceMap *outmap = NULL;
250 double distance = 0;
251 class BoundaryTriangleSet *triangle = NULL;
252 Vector centroid;
253
254 if ((Surface == NULL) || (LC == NULL) || (molecules->ListOfMolecules.empty())) {
255 DoeLog(1) && (eLog()<< Verbose(1) <<"No Tesselation, no LinkedCell or no molecule given." << endl);
256 return outmap;
257 }
258 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
259 (*MolWalker)->doCountAtoms();
260 outmap = new CorrelationToSurfaceMap;
261 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
262 if ((*MolWalker)->ActiveFlag) {
263 DoLog(1) && (Log() << Verbose(1) << "Current molecule is " << (*MolWalker)->name << "." << endl);
264 if ((*MolWalker)->empty())
265 DoLog(1) && (1) && (Log() << Verbose(1) << "\t is empty." << endl);
266 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
267 DoLog(1) && (Log() << Verbose(1) << "\tCurrent atom is " << *(*iter) << "." << endl);
268 if ((type == NULL) || ((*iter)->type == type)) {
269 TriangleIntersectionList Intersections((*iter)->node,Surface,LC);
270 distance = Intersections.GetSmallestDistance();
271 triangle = Intersections.GetClosestTriangle();
272 outmap->insert ( pair<double, pair<atom *, BoundaryTriangleSet*> >(distance, pair<atom *, BoundaryTriangleSet*> ((*iter), triangle) ) );
273 }
274 }
275 } else {
276 DoLog(1) && (Log() << Verbose(1) << "molecule " << (*MolWalker)->name << " is not active." << endl);
277 }
278
279 return outmap;
280};
281
282/** Calculates the distance (pair) correlation between a given element, all its periodic images and and a surface.
283 * Note that we also put all periodic images found in the cells given by [ -ranges[i], ranges[i] ] and i=0,...,NDIM-1.
284 * I.e. We multiply the atom::node with the inverse of the domain matrix, i.e. transform it to \f$[0,0^3\f$, then add per
285 * axis an integer from [ -ranges[i], ranges[i] ] onto it and multiply with the domain matrix to bring it back into
286 * the real space. Then, we Tesselation::FindClosestTriangleToPoint() and DistanceToTrianglePlane().
287 * \param *out output stream for debugging
288 * \param *molecules list of molecules structure
289 * \param *type element or NULL (if any element)
290 * \param *Surface pointer to Tesselation class surface
291 * \param *LC LinkedCell structure to quickly find neighbouring atoms
292 * \param ranges[NDIM] interval boundaries for the periodic images to scan also
293 * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest
294 */
295CorrelationToSurfaceMap *PeriodicCorrelationToSurface(MoleculeListClass * const &molecules, const element * const type, const Tesselation * const Surface, const LinkedCell *LC, const int ranges[NDIM] )
296{
297 Info FunctionInfo(__func__);
298 CorrelationToSurfaceMap *outmap = NULL;
299 double distance = 0;
300 class BoundaryTriangleSet *triangle = NULL;
301 Vector centroid;
302 int n[NDIM];
303 Vector periodicX;
304 Vector checkX;
305
306 if ((Surface == NULL) || (LC == NULL) || (molecules->ListOfMolecules.empty())) {
307 DoLog(1) && (Log() << Verbose(1) <<"No Tesselation, no LinkedCell or no molecule given." << endl);
308 return outmap;
309 }
310 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
311 (*MolWalker)->doCountAtoms();
312 outmap = new CorrelationToSurfaceMap;
313 double ShortestDistance = 0.;
314 BoundaryTriangleSet *ShortestTriangle = NULL;
315 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin(); MolWalker != molecules->ListOfMolecules.end(); MolWalker++)
316 if ((*MolWalker)->ActiveFlag) {
317 double * FullMatrix = ReturnFullMatrixforSymmetric(World::getInstance().getDomain());
318 double * FullInverseMatrix = InverseMatrix(FullMatrix);
319 DoLog(2) && (Log() << Verbose(2) << "Current molecule is " << *MolWalker << "." << endl);
320 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
321 DoLog(3) && (Log() << Verbose(3) << "Current atom is " << **iter << "." << endl);
322 if ((type == NULL) || ((*iter)->type == type)) {
323 periodicX = *(*iter)->node;
324 periodicX.MatrixMultiplication(FullInverseMatrix); // x now in [0,1)^3
325 // go through every range in xyz and get distance
326 ShortestDistance = -1.;
327 for (n[0]=-ranges[0]; n[0] <= ranges[0]; n[0]++)
328 for (n[1]=-ranges[1]; n[1] <= ranges[1]; n[1]++)
329 for (n[2]=-ranges[2]; n[2] <= ranges[2]; n[2]++) {
330 checkX = Vector(n[0], n[1], n[2]) + periodicX;
331 checkX.MatrixMultiplication(FullMatrix);
332 TriangleIntersectionList Intersections(&checkX,Surface,LC);
333 distance = Intersections.GetSmallestDistance();
334 triangle = Intersections.GetClosestTriangle();
335 if ((ShortestDistance == -1.) || (distance < ShortestDistance)) {
336 ShortestDistance = distance;
337 ShortestTriangle = triangle;
338 }
339 }
340 // insert
341 outmap->insert ( pair<double, pair<atom *, BoundaryTriangleSet*> >(ShortestDistance, pair<atom *, BoundaryTriangleSet*> (*iter, ShortestTriangle) ) );
342 //Log() << Verbose(1) << "INFO: Inserting " << Walker << " with distance " << ShortestDistance << " to " << *ShortestTriangle << "." << endl;
343 }
344 }
345 delete[](FullMatrix);
346 delete[](FullInverseMatrix);
347 }
348
349 return outmap;
350};
351
352/** Returns the index of the bin for a given value.
353 * \param value value whose bin to look for
354 * \param BinWidth width of bin
355 * \param BinStart first bin
356 */
357int GetBin ( const double value, const double BinWidth, const double BinStart )
358{
359 Info FunctionInfo(__func__);
360 int bin =(int) (floor((value - BinStart)/BinWidth));
361 return (bin);
362};
363
364
365/** Prints correlation (double, int) pairs to file.
366 * \param *file file to write to
367 * \param *map map to write
368 */
369void OutputCorrelation( ofstream * const file, const BinPairMap * const map )
370{
371 Info FunctionInfo(__func__);
372 *file << "BinStart\tCount" << endl;
373 for (BinPairMap::const_iterator runner = map->begin(); runner != map->end(); ++runner) {
374 *file << setprecision(8) << runner->first << "\t" << runner->second << endl;
375 }
376};
377
378/** Prints correlation (double, (atom*,atom*) ) pairs to file.
379 * \param *file file to write to
380 * \param *map map to write
381 */
382void OutputPairCorrelation( ofstream * const file, const PairCorrelationMap * const map )
383{
384 Info FunctionInfo(__func__);
385 *file << "BinStart\tAtom1\tAtom2" << endl;
386 for (PairCorrelationMap::const_iterator runner = map->begin(); runner != map->end(); ++runner) {
387 *file << setprecision(8) << runner->first << "\t" << *(runner->second.first) << "\t" << *(runner->second.second) << endl;
388 }
389};
390
391/** Prints correlation (double, int) pairs to file.
392 * \param *file file to write to
393 * \param *map map to write
394 */
395void OutputCorrelationToPoint( ofstream * const file, const CorrelationToPointMap * const map )
396{
397 Info FunctionInfo(__func__);
398 *file << "BinStart\tAtom::x[i]-point.x[i]" << endl;
399 for (CorrelationToPointMap::const_iterator runner = map->begin(); runner != map->end(); ++runner) {
400 *file << runner->first;
401 for (int i=0;i<NDIM;i++)
402 *file << "\t" << setprecision(8) << (runner->second.first->node->at(i) - runner->second.second->at(i));
403 *file << endl;
404 }
405};
406
407/** Prints correlation (double, int) pairs to file.
408 * \param *file file to write to
409 * \param *map map to write
410 */
411void OutputCorrelationToSurface( ofstream * const file, const CorrelationToSurfaceMap * const map )
412{
413 Info FunctionInfo(__func__);
414 *file << "BinStart\tTriangle" << endl;
415 if (!map->empty())
416 for (CorrelationToSurfaceMap::const_iterator runner = map->begin(); runner != map->end(); ++runner) {
417 *file << setprecision(8) << runner->first << "\t" << *(runner->second.first) << "\t" << *(runner->second.second) << endl;
418 }
419};
420
Note: See TracBrowser for help on using the repository browser.