generic.hpp
Go to the documentation of this file.
1 //-----------------------------------------------------------------------------------
2 // d-SEAMS - Deferred Structural Elucidation Analysis for Molecular Simulations
3 //
4 // Copyright (c) 2018--present d-SEAMS core team
5 //
6 // This program is free software: you can redistribute it and/or modify
8 // the Open Source Initiative.
9 //
10 // A copy of the MIT License is included in the LICENSE file of this repository.
11 // You should have received a copy of the MIT License along with this program.
12 // If not, see <https://opensource.org/licenses/MIT>.
13 //-----------------------------------------------------------------------------------
14
15 #ifndef __GENERIC_H_
16 #define __GENERIC_H_
17
18 #include <array>
19 #include <iostream>
20 #include <math.h>
21 #include <mol_sys.hpp>
22
23 // Boost
24 #include <boost/math/constants/constants.hpp>
25 // Eigen
26 #include <Eigen/Core>
27 #include <Eigen/Dense>
28
49 namespace gen {
50
54 const double pi = boost::math::constants::pi<double>();
55
61 inline double radDeg(double angle) { return (angle * 180) / gen::pi; }
62
65
68
77 inline double calcMedian(std::vector<double> *input) {
78  int n = (*input).size(); // Number of elements
79  double median; // Output median value
80
81  // Sort the vector
82  std::sort((*input).begin(), (*input).end());
83
84  // Calculate the median
85  // For even values, the median is the average of the two middle values
86  if (n % 2 == 0) {
87  median = 0.5 * ((*input)[n / 2] + (*input)[n / 2 - 1]); // n/2+n/2-1
88  } // median is average of middle values
89  else {
90  median = (*input)[(n + 1) / 2 -
91  1]; // middle value of 7 elements is the 4th element
92  } // if odd, it is the middle value
93
94  return median;
95 }
96
97 // Generic function for getting the unwrapped distance
107 inline double
109  int iatom, int jatom) {
111  double r2 = 0.0; // Squared absolute distance
112
113  // Get x1-x2 etc
114  dr[0] = fabs(yCloud->pts[iatom].x - yCloud->pts[jatom].x);
115  dr[1] = fabs(yCloud->pts[iatom].y - yCloud->pts[jatom].y);
116  dr[2] = fabs(yCloud->pts[iatom].z - yCloud->pts[jatom].z);
117
118  // Get the squared absolute distance
119  for (int k = 0; k < 3; k++) {
120  // Correct for periodicity
121  dr[k] -= yCloud->box[k] * round(dr[k] / yCloud->box[k]);
122  r2 += pow(dr[k], 2.0);
123  }
124
125  return sqrt(r2);
126 }
127
139  molSys::PointCloud<molSys::Point<double>, double> *yCloud, int iatom,
140  std::vector<double> singlePoint) {
142  double r2 = 0.0; // Squared absolute distance
143
144  // Get x1-x2 etc
145  dr[0] = fabs(yCloud->pts[iatom].x - singlePoint[0]);
146  dr[1] = fabs(yCloud->pts[iatom].y - singlePoint[1]);
147  dr[2] = fabs(yCloud->pts[iatom].z - singlePoint[2]);
148
149  // Get the squared absolute distance
150  for (int k = 0; k < 3; k++) {
151  // Correct for periodicity
152  dr[k] -= yCloud->box[k] * round(dr[k] / yCloud->box[k]);
153  r2 += pow(dr[k], 2.0);
154  }
155
156  return sqrt(r2);
157 }
158
159 // Generic function for getting the distance (no PBCs applied)
169 inline double
171  int jatom) {
173  double r2 = 0.0; // Squared absolute distance
174
175  // Get x1-x2 etc
176  dr[0] = fabs(yCloud->pts[iatom].x - yCloud->pts[jatom].x);
177  dr[1] = fabs(yCloud->pts[iatom].y - yCloud->pts[jatom].y);
178  dr[2] = fabs(yCloud->pts[iatom].z - yCloud->pts[jatom].z);
179
180  // Get the squared absolute distance
181  for (int k = 0; k < 3; k++) {
182  r2 += pow(dr[k], 2.0);
183  }
184
185  return sqrt(r2);
186 }
187
188 // Generic function for getting the relative coordinates
201  int jatom) {
203  std::array<double, 3> box = {yCloud->box[0], yCloud->box[1], yCloud->box[2]};
204  double r2 = 0.0; // Squared absolute distance
205
206  // Get x1-x2 etc
207  dr[0] = yCloud->pts[iatom].x - yCloud->pts[jatom].x;
208  dr[1] = yCloud->pts[iatom].y - yCloud->pts[jatom].y;
209  dr[2] = yCloud->pts[iatom].z - yCloud->pts[jatom].z;
210
211  // Get the relative distance
212  for (int k = 0; k < 3; k++) {
213  //
214  if (dr[k] < -box[k] * 0.5) {
215  dr[k] = dr[k] + box[k];
216  }
217  if (dr[k] >= box[k] * 0.5) {
218  dr[k] = dr[k] - box[k];
219  }
220  }
221
222  return dr;
223 }
224
225 // Function for sorting according to atom ID
226 // Comparator for std::sort
235  const molSys::Point<double> &b) {
236  return a.atomID < b.atomID;
237 }
238
241  std::string outFile);
242
245  molSys::PointCloud<molSys::Point<double>, double> *yCloud, int iatomIndex,
246  int jatomIndex, double *x_i, double *y_i, double *z_i, double *x_j,
247  double *y_j, double *z_j);
248
252  std::vector<double> quat2);
253
260  std::istringstream iss(line);
264 }
265
271  std::istringstream iss(line);
272  std::vector<double> tokens;
273  double number; // Each number being read in from the line
274  while (iss >> number) {
275  tokens.push_back(number);
276  }
278 }
279
285  std::istringstream iss(line);
286  std::vector<int> tokens;
287  int number; // Each number being read in from the line
288  while (iss >> number) {
289  tokens.push_back(number);
290  }
292 }
293
298 inline bool file_exists(const std::string &name) {
299  // Replace by boost function later
300  struct stat buffer;
301  return (stat(name.c_str(), &buffer) == 0);
302 }
303
314  if (neigh == 0) {
315  return v;
316  }
317  for (int m = 0; m < 2 * l + 1; m++) {
318  v[m] = (1.0 / (double)neigh) * v[m];
319  }
320
321  return v;
322 }
323
324 } // namespace gen
325
326 #endif // __NEIGHBOURS_H_
T c_str(T... args)
double distance(molSys::PointCloud< molSys::Point< double >, double > *yCloud, int iatom, int jatom)
Inline generic function for obtaining the wrapped distance between two particles WITHOUT applying PBC...
Definition: generic.hpp:170
std::vector< std::complex< double > > avgVector(std::vector< std::complex< double >> v, int l, int neigh)
Definition: generic.hpp:313
std::array< double, 3 > relDist(molSys::PointCloud< molSys::Point< double >, double > *yCloud, int iatom, int jatom)
Definition: generic.hpp:200
double eigenVecAngle(std::vector< double > OO, std::vector< double > OH)
Eigen function for getting the angle (in radians) between the O–O and O-H vectors.
Definition: generic.cpp:155
double calcMedian(std::vector< double > *input)
Inline generic function for calculating the median given a vector of the values.
Definition: generic.hpp:77
std::vector< int > tokenizerInt(std::string line)
Function for tokenizing line strings into a vector of ints.
Definition: generic.hpp:284
double periodicDist(molSys::PointCloud< molSys::Point< double >, double > *yCloud, int iatom, int jatom)
Inline generic function for obtaining the unwrapped periodic distance between two particles,...
Definition: generic.hpp:108
double unWrappedDistFromPoint(molSys::PointCloud< molSys::Point< double >, double > *yCloud, int iatom, std::vector< double > singlePoint)
Definition: generic.hpp:138
bool file_exists(const std::string &name)
Function for checking if a file exists or not.
Definition: generic.hpp:298
const double pi
Definition: generic.hpp:54
std::vector< std::string > tokenizer(std::string line)
Function for tokenizing line strings into words (strings) delimited by whitespace....
Definition: generic.hpp:259
double angDistDegQuaternions(std::vector< double > quat1, std::vector< double > quat2)
Definition: generic.cpp:267
int unwrappedCoordShift(molSys::PointCloud< molSys::Point< double >, double > *yCloud, int iatomIndex, int jatomIndex, double *x_i, double *y_i, double *z_i, double *x_j, double *y_j, double *z_j)
Shift particles (unwrapped coordinates)
Definition: generic.cpp:71
Definition: generic.hpp:61
bool compareByAtomID(const molSys::Point< double > &a, const molSys::Point< double > &b)
Definition: generic.hpp:234
int prettyPrintYoda(molSys::PointCloud< molSys::Point< double >, double > *yCloud, std::string outFile)
Generic function for printing all the struct information.
Definition: generic.cpp:25
double getAverageWithoutOutliers(std::vector< double > inpVec)
Get the average, after excluding the outliers, using quartiles.
Definition: generic.cpp:169
std::vector< double > tokenizerDouble(std::string line)
Function for tokenizing line strings into a vector of doubles.
Definition: generic.hpp:270
The main molecular system handler.
Small generic functions that are shared by all namespaces.
Definition: generic.hpp:49
T push_back(T... args)
T sort(T... args)
This contains a collection of points; contains information for a particular frame.
Definition: mol_sys.hpp:170
This contains per-particle information.
Definition: mol_sys.hpp:149