ProSHADE  0.6.5 (NOV 2018)
Protein Shape Descriptors and Symmetry Detection
ProSHADE_descriptors.cpp
Go to the documentation of this file.
1 
21 //============================================ FFTW3 + SOFT
22 #ifdef __cplusplus
23 extern "C" {
24 #endif
25 #include <fftw3.h>
26 #include <wrap_fftw.h>
27 #include <makeweights.h>
28 #include <s2_primitive.h>
29 #include <s2_cospmls.h>
30 #include <s2_legendreTransforms.h>
31 #include <s2_semi_fly.h>
32 #include <rotate_so3_utils.h>
33 #include <utils_so3.h>
34 #include <soft_fftw.h>
35 #include <rotate_so3_fftw.h>
36 #ifdef __cplusplus
37 }
38 #endif
39 
40 //============================================ RVAPI
41 #include <rvapi_interface.h>
42 
43 //============================================ ProSHADE
44 #include "ProSHADE.h"
45 #include "ProSHADE_internal.h"
46 #include "ProSHADE_misc.h"
47 #include "ProSHADE_legendre.h"
48 
62 double ProSHADE_internal_misc::pearsonCorrCoeff ( double *valSet1,
63  double *valSet2,
64  unsigned int length )
65 {
66  //======================================== Find vector means
67  double xMean = 0.0;
68  double yMean = 0.0;
69  for ( unsigned int iter = 0; iter < length; iter++ )
70  {
71  xMean += valSet1[iter];
72  yMean += valSet2[iter];
73  }
74  xMean /= static_cast<double> ( length );
75  yMean /= static_cast<double> ( length );
76 
77  //======================================== Get Pearson's correlation coefficient
78  double xmmymm = 0.0;
79  double xmmsq = 0.0;
80  double ymmsq = 0.0;
81  for ( unsigned int iter = 0; iter < length; iter++ )
82  {
83  xmmymm += ( valSet1[iter] - xMean ) * ( valSet2[iter] - yMean );
84  xmmsq += pow( valSet1[iter] - xMean, 2.0 );
85  ymmsq += pow( valSet2[iter] - yMean, 2.0 );
86  }
87 
88  //======================================== Done
89  double ret = xmmymm / ( sqrt(xmmsq) * sqrt(ymmsq) );
90  if ( std::isnan (ret) )
91  {
92  return ( 0.0 );
93  }
94  return ( ret );
95 
96 }
97 
112 std::array<double,2> ProSHADE_internal_misc::complexMultiplication ( double* r1,
113  double* i1,
114  double* r2,
115  double* i2 )
116 {
117  //======================================== Initialisation
118  std::array<double,2> ret;
119 
120  //======================================== Multiplication
121  ret[0] = (*r1)*(*r2) - (*i1)*(*i2);
122  ret[1] = (*r1)*(*i2) + (*i1)*(*r2);
123 
124  //======================================== Return
125  return ( ret );
126 
127 }
128 
143 std::array<double,2> ProSHADE_internal_misc::complexMultiplicationConjug ( double* r1,
144  double* i1,
145  double* r2,
146  double* i2 )
147 {
148  //======================================== Initialisation
149  std::array<double,2> ret;
150 
151  //======================================== Multiplication
152  ret[0] = (*r1)*(*r2) + (*i1)*(*i2);
153  ret[1] = -(*r1)*(*i2) + (*i1)*(*r2);
154 
155  //======================================== Return
156  return ( ret );
157 
158 }
159 
174 double ProSHADE_internal_misc::complexMultiplicationReal ( double* r1,
175  double* i1,
176  double* r2,
177  double* i2 )
178 {
179  //======================================== Return
180  return ( (*r1)*(*r2) - (*i1)*(*i2) );
181 
182 }
183 
199 double ProSHADE_internal_misc::complexMultiplicationConjugReal ( double* r1,
200  double* i1,
201  double* r2,
202  double* i2 )
203 {
204  //======================================== Return
205  return ( (*r1)*(*r2) + (*i1)*(*i2) );
206 
207 }
208 
219 void ProSHADE_internal::ProSHADE_data::precomputeRotInvDescriptor ( ProSHADE::ProSHADE_settings* settings )
220 {
221  //======================================== Sanity check
222  if ( !this->_sphericalCoefficientsComputed )
223  {
224  std::cerr << "!!! ProSHADE ERROR !!! Error in file " << this->_inputFileName << " !!! Attempted to pre-compute RRP matrices for RotInv descriptor before computing the spherical harmonics decomposition. Call the getSphericalHarmonicsCoeffs function BEFORE the precomputeRotInvDescriptor function." << std::endl;
225 
226  if ( settings->htmlReport )
227  {
228  std::stringstream hlpSS;
229  hlpSS << "<font color=\"red\">" << "Attempted to pre-compute RRP matrices for RotInv descriptor before computing the spherical harmonics decomposition. This looks like an internal bug, please report this case." << "</font>";
230  rvapi_set_text ( hlpSS.str().c_str(),
231  "ProgressSection",
232  settings->htmlReportLineProgress,
233  1,
234  1,
235  1 );
236  settings->htmlReportLineProgress += 1;
237  rvapi_flush ( );
238  }
239 
240  exit ( -1 );
241  }
242 
243  //======================================== Initialise internal values
244  this->_rrpMatrices = new double** [this->_bandwidthLimit];
245  for ( unsigned int bwIt = 0; bwIt < this->_bandwidthLimit; bwIt++ )
246  {
247  //==================================== Odd bands are 0, so just ignore them
248  if ( !this->_keepOrRemove ) { if ( ( bwIt % 2 ) != 0 ) { continue; } }
249 
250  this->_rrpMatrices[bwIt] = new double* [this->_noShellsWithData];
251  for ( unsigned int shIt = 0; shIt < this->_noShellsWithData; shIt++ )
252  {
253  this->_rrpMatrices[bwIt][shIt] = new double [this->_noShellsWithData];
254  }
255  }
256 
257  //======================================== Initialise local variables
258  double descValR = 0.0;
259  unsigned int arrPos = 0;
260 
261  //======================================== For each band (l)
262  for ( unsigned int band = 0; band < this->_bandwidthLimit; band++ )
263  {
264  //==================================== Ignore odd bands, they are always 0
265  if ( !this->_keepOrRemove ) { if ( band % 2 == 1 ) { continue; } }
266 
267  //==================================== For each unique shell couple
268  for ( unsigned int shell1 = 0; shell1 < this->_noShellsWithData; shell1++ )
269  {
270  for ( unsigned int shell2 = 0; shell2 < this->_noShellsWithData; shell2++ )
271  {
272  //============================ Make use of symmetry
273  if ( shell1 > shell2 ) { continue; }
274 
275  //============================ Initialise
276  descValR = 0.0;
277 
278  //============================ Sum over order (m)
279  for ( unsigned int order = 0; order < static_cast<unsigned int> ( ( 2 * band ) + 1 ); order++ )
280  {
281  arrPos = seanindex ( order-band, band, this->_bandwidthLimit );
282  descValR += ( this->_realSHCoeffs[shell1][arrPos] * this->_realSHCoeffs[shell2][arrPos] ) +
283  ( this->_imagSHCoeffs[shell1][arrPos] * this->_imagSHCoeffs[shell2][arrPos] );
284  }
285 
286  //============================ Save the matrices
287  this->_rrpMatrices[band][shell1][shell2] = descValR;
288  this->_rrpMatrices[band][shell2][shell1] = descValR;
289  }
290  }
291  }
292 
293  //======================================== Done
294  this->_rrpMatricesPrecomputed = true;
295  return ;
296 
297 }
298 
310 std::vector<double> ProSHADE_internal::ProSHADE_compareOneAgainstAll::getEnergyLevelsDistance ( int verbose, ProSHADE::ProSHADE_settings* settings )
311 {
312  //======================================== Sanity check
313  if ( !this->one->_rrpMatricesPrecomputed )
314  {
315  std::cerr << "!!! ProSHADE ERROR !!! Error in file compatison !!! Attempted to get energy level distance for structure objects where at least one of them did not have the RotInv pre-computed. Call the precomputeRotInvDescriptor on all structure objects before giving them to ProSHADE_compareOneAgainstAll constructor." << std::endl;
316 
317  if ( settings->htmlReport )
318  {
319  std::stringstream hlpSS;
320  hlpSS << "<font color=\"red\">" << "Attempted to pre-compute RRP matrices for RotInv descriptor before computing the spherical harmonics decomposition. This looks like an internal bug, please report this case." << "</font>";
321  rvapi_set_text ( hlpSS.str().c_str(),
322  "ProgressSection",
323  settings->htmlReportLineProgress,
324  1,
325  1,
326  1 );
327  settings->htmlReportLineProgress += 1;
328  rvapi_flush ( );
329  }
330 
331  exit ( -1 );
332  }
333  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( this->all->size() ); iter++ )
334  {
335  if ( !this->all->at(iter)->_rrpMatricesPrecomputed )
336  {
337  std::cerr << "!!! ProSHADE ERROR !!! Error in file compatison !!! Attempted to get energy level distance for structure objects where at least one of them did not have the RotInv pre-computed. Call the precomputeRotInvDescriptor on all structure objects before giving them to ProSHADE_compareOneAgainstAll constructor." << std::endl;
338 
339  if ( settings->htmlReport )
340  {
341  std::stringstream hlpSS;
342  hlpSS << "<font color=\"red\">" << "Attempted to get energy level distance for structure objects where at least one of them did not have the RotInv pre-computed. This looks like an internal bug, please report this case." << "</font>";
343  rvapi_set_text ( hlpSS.str().c_str(),
344  "ProgressSection",
345  settings->htmlReportLineProgress,
346  1,
347  1,
348  1 );
349  settings->htmlReportLineProgress += 1;
350  rvapi_flush ( );
351  }
352 
353  exit ( -1 );
354  }
355  }
356  if ( verbose > 3 )
357  {
358  std::cout << ">>>>>>>> Sanity checks passed." << std::endl;
359  }
360 
361  //======================================== If already done, do not repeat :-)
362  if ( this->_energyLevelsComputed )
363  {
364  return ( this->_distancesEnergyLevels );
365  }
366 
367  //======================================== Initialise local variables
368  unsigned int arrIter = 0;
369  bool ignoreThisL = false;
370  std::vector<double> bandDists;
371  this->_distancesEnergyLevels = std::vector<double> ( static_cast<int> ( this->all->size() ) );
372 
373  for ( unsigned int strIt = 0; strIt < static_cast<unsigned int> ( this->all->size() ); strIt++ )
374  {
375  if ( verbose > 3 )
376  {
377  std::cout << ">>>>>>>> Now computing distance for structure 0 against structure " << strIt+1 << " ." << std::endl;
378  }
379 
380  //==================================== Clean
381  bandDists.clear ( );
382 
383  //==================================== Set up
384  unsigned int minShellsToUse = std::min ( this->one->_noShellsWithData, this->all->at(strIt)->_noShellsWithData );
385  unsigned int bandLimit = std::min ( this->one->_bandwidthLimit, this->all->at(strIt)->_bandwidthLimit );
386 
387  //==================================== Initialise local variables
388  double *str1Vals = new double[minShellsToUse * minShellsToUse];
389  double *str2Vals = new double[minShellsToUse * minShellsToUse];
390 
391  //==================================== For each band (l)
392  for ( unsigned int band = 0; band < bandLimit; band++ )
393  {
394  //================================ Check if this band is to be used
395  if ( !this->_keepOrRemove ) { if ( band % 2 != 0 ) { continue; } }
396 
397  ignoreThisL = false;
398  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) == band ) { ignoreThisL = true; } }
399  if ( ignoreThisL ) { continue; }
400 
401  //================================ Reset local counter
402  arrIter = 0;
403 
404  //================================ For each shell pair
405  for ( unsigned int shell1 = 0; shell1 < minShellsToUse; shell1++ )
406  {
407  for ( unsigned int shell2 = 0; shell2 < minShellsToUse; shell2++ )
408  {
409  str1Vals[arrIter] = this->one->_rrpMatrices[band][shell1][shell2] *
410  pow ( static_cast<double> ( shell1 ), this->_matrixPowerWeight ) *
411  pow ( static_cast<double> ( shell2 ), this->_matrixPowerWeight );
412  str2Vals[arrIter] = this->all->at(strIt)->_rrpMatrices[band][shell1][shell2] *
413  pow ( static_cast<double> ( shell1 ), this->_matrixPowerWeight ) *
414  pow ( static_cast<double> ( shell2 ), this->_matrixPowerWeight );
415 
416  arrIter += 1;
417  }
418  }
419 
420  //================================ Get Pearson's Correlation Coefficient
421  bandDists.emplace_back ( ProSHADE_internal_misc::pearsonCorrCoeff ( str1Vals, str2Vals, ( minShellsToUse * minShellsToUse ) ) );
422  }
423 
424  //==================================== Get distance
425  this->_distancesEnergyLevels.at( strIt ) = static_cast<double> ( std::accumulate ( bandDists.begin(), bandDists.end(), 0.0 ) ) / static_cast<double> ( bandDists.size() );
426 
427  //==================================== Clean up
428  delete[] str1Vals;
429  delete[] str2Vals;
430 
431  if ( verbose > 2 )
432  {
433  std::cout << ">>>>> Distance for structure 0 against structure " << strIt+1 << " complete." << std::endl;
434  }
435  }
436 
437  //======================================== Done
438  this->_energyLevelsComputed = true;
439 
440  //======================================== Return
441  return ( this->_distancesEnergyLevels );
442 
443 }
444 
454 void ProSHADE_internal::ProSHADE_comparePairwise::precomputeTrSigmaDescriptor ( )
455 {
456  //======================================== Figure required number of bands
457  unsigned int numberBandsToUse = 0;
458  bool toBeIgnored = false;
459 
460  for ( unsigned int bandIter = 0; bandIter < this->_bandwidthLimit; bandIter++ )
461  {
462  //==================================== Ignore odd l's as they are 0 anyway
463  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
464 
465  //==================================== Ignore l's designated to be ignored
466  toBeIgnored = false;
467  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
468  if ( toBeIgnored ) { continue; }
469 
470  //==================================== This band will be used!
471  numberBandsToUse += 1;
472  }
473 
474  //======================================== Initialise local variables
475  std::vector<double> radiiVals1;
476  std::vector<double> radiiVals2;
477  std::vector< std::array<double,2> > radiiVals;
478  std::array<double,2> arrVal1;
479  unsigned int index1 = 0;
480  unsigned int index2 = 0;
481  unsigned int bandUsageIndex = 0;
482  double correlNormFactor = 0.0;
483 
484  //======================================== Initialise internal values
485  this->_trSigmaWeights = std::vector<double> ( 2 );
486  this->_trSigmaWeights.at(0) = 0.0;
487  this->_trSigmaWeights.at(1) = 0.0;
488 
489  this->_trSigmaEMatrix = new std::array<double,2>** [numberBandsToUse];
490  for ( unsigned int bandIter = 0; bandIter < this->_bandwidthLimit; bandIter++ )
491  {
492  //==================================== Ignore odd l's as they are 0 anyway
493  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
494 
495  //==================================== Ignore l's designated to be ignored
496  toBeIgnored = false;
497  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
498  if ( toBeIgnored ) { continue; }
499 
500  this->_trSigmaEMatrix[bandUsageIndex] = new std::array<double,2>* [static_cast<unsigned int> ( ( bandIter * 2 ) + 1 )];
501  for ( unsigned int lIter = 0; lIter < static_cast<unsigned int> ( ( bandIter * 2 ) + 1 ); lIter++ )
502  {
503  this->_trSigmaEMatrix[bandUsageIndex][lIter] = new std::array<double,2> [static_cast<unsigned int> ( ( bandIter * 2 ) + 1 )];
504  }
505 
506  bandUsageIndex += 1;
507  }
508 
509  //======================================== For each band (l)
510  bandUsageIndex = 0;
511  for ( unsigned int bandIter = 0; bandIter < this->_bandwidthLimit; bandIter++ )
512  {
513  //==================================== Ignore odd l's as they are 0 anyway
514  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
515 
516  //==================================== Ignore l's designated to be ignored
517  toBeIgnored = false;
518  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
519  if ( toBeIgnored ) { continue; }
520 
521  //==================================== For each order (m)
522  for ( unsigned int orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
523  {
524  //================================ Clear values for weights
525  radiiVals1.clear ( );
526  radiiVals2.clear ( );
527  radiiVals1 = std::vector<double> ( this->_noShellsObj1 );
528  radiiVals2 = std::vector<double> ( this->_noShellsObj2 );
529 
530  //================================ Get index for m
531  index1 = seanindex ( orderIter-bandIter, bandIter, this->_bandwidthLimit );
532 
533  //================================ For each radius, deal with weights
534  for ( unsigned int radiusIter = 0; radiusIter < this->_maxShellsToUse; radiusIter++ )
535  {
536  if ( radiusIter < this->_noShellsObj1 )
537  {
538  //======================== Multiply with itself - i.e. square
539  radiiVals1.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjugReal ( &this->_obj1RealCoeffs[radiusIter][index1],
540  &this->_obj1ImagCoeffs[radiusIter][index1],
541  &this->_obj1RealCoeffs[radiusIter][index1],
542  &this->_obj1ImagCoeffs[radiusIter][index1] );
543 
544  //======================== Apply the integral r^2
545  radiiVals1.at(radiusIter)*= pow ( ( static_cast<double> ( radiusIter + 1 ) * this->_shellSpacing ), 2.0 );
546  }
547 
548  if ( radiusIter < this->_noShellsObj2 )
549  {
550  //======================== Multiply with itself - i.e. square
551  radiiVals2.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjugReal ( &this->_obj2RealCoeffs[radiusIter][index1],
552  &this->_obj2ImagCoeffs[radiusIter][index1],
553  &this->_obj2RealCoeffs[radiusIter][index1],
554  &this->_obj2ImagCoeffs[radiusIter][index1] );
555 
556  //======================== Apply the integral r^2
557  radiiVals2.at(radiusIter)*= pow ( ( static_cast<double> ( radiusIter + 1 ) * this->_shellSpacing ), 2.0 );
558  }
559  }
560 
561  //================================ Integrate weights
562  this->_trSigmaWeights.at(0) += this->gl20IntRR ( &radiiVals1 );
563  this->_trSigmaWeights.at(1) += this->gl20IntRR ( &radiiVals2 );
564 
565  //================================ For each combination of m and m' for E matrices
566  for ( unsigned int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
567  {
568  //============================ Set index for order 2 (i.e. m' in E_{l,m,m'})
569  index2 = seanindex ( order2Iter-bandIter, bandIter, this->_bandwidthLimit );
570 
571  //============================ Clear the vector to integrate over
572  radiiVals.clear ( );
573  radiiVals = std::vector< std::array<double,2> > ( this->_minShellsToUse );
574 
575  //============================ Find the c*c values for different radii
576  for ( unsigned int radiusIter = 0; radiusIter < this->_minShellsToUse; radiusIter++ )
577  {
578  //======================== Multiply coeffs
579  radiiVals.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjug ( &this->_obj1RealCoeffs[radiusIter][index1],
580  &this->_obj1ImagCoeffs[radiusIter][index1],
581  &this->_obj2RealCoeffs[radiusIter][index2],
582  &this->_obj2ImagCoeffs[radiusIter][index2] );
583 
584  //======================== Apply r^2 integral weight
585  radiiVals.at(radiusIter)[0] *= pow ( ( static_cast<double> ( radiusIter + 1 ) * this->_shellSpacing ), 2.0 );
586  radiiVals.at(radiusIter)[1] *= pow ( ( static_cast<double> ( radiusIter + 1 ) * this->_shellSpacing ), 2.0 );
587  }
588 
589  //============================ Integrate over all radii using X-point Gauss-Legendre
590  arrVal1 = this->gl20IntCR ( &radiiVals );
591 
592  //============================ Save the result into E matrices
593  this->_trSigmaEMatrix[bandUsageIndex][orderIter][order2Iter][0] = arrVal1[0];
594  this->_trSigmaEMatrix[bandUsageIndex][orderIter][order2Iter][1] = arrVal1[1];
595  }
596  }
597  bandUsageIndex += 1;
598  }
599 
600  //======================================== For each band (l), normalise using the correlation-like formula
601  bandUsageIndex = 0;
602  correlNormFactor = sqrt ( this->_trSigmaWeights.at(0) * this->_trSigmaWeights.at(1) );
603 
604  for ( unsigned int bandIter = 0; bandIter < this->_bandwidthLimit; bandIter++ )
605  {
606  //==================================== Ignore odd l's as they are 0 anyway
607  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
608 
609  //==================================== Ignore l's designated to be ignored
610  toBeIgnored = false;
611  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
612  if ( toBeIgnored ) { continue; }
613 
614  //==================================== For each combination of m and m' for E matrices
615  for ( unsigned int order1Iter = 0; order1Iter < ( ( bandIter * 2 ) + 1 ); order1Iter++ )
616  {
617  for ( unsigned int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
618  {
619  //============================ Normalisation using the correlation formula
620  this->_trSigmaEMatrix[bandUsageIndex][order1Iter][order2Iter][0] /= correlNormFactor;
621  this->_trSigmaEMatrix[bandUsageIndex][order1Iter][order2Iter][1] /= correlNormFactor;
622  }
623  }
624 
625  bandUsageIndex += 1;
626  }
627 
628  //======================================== Done
629  this->_trSigmaPreComputed = true;
630 
631 }
632 
646 void ProSHADE_internal::ProSHADE_compareOneAgainstAll::precomputeTrSigmaDescriptor ( double shellSpacing, std::vector<unsigned int>* glIntegOrderVec, ProSHADE::ProSHADE_settings* settings )
647 {
648  //======================================== Initialise local variables
649  std::vector<double> radiiVals1;
650  std::vector<double> radiiVals2;
651  std::vector< std::array<double,2> > radiiVals;
652  std::array<double,2> arrVal1;
653  std::vector<double> glAbscissas;
654  std::vector<double> glWeights;
655  unsigned int glIntegrationOrder = 0;
656  unsigned int index11 = 0;
657  unsigned int index12 = 0;
658  unsigned int index21 = 0;
659  unsigned int index22 = 0;
660  double correlNormFactor = 0.0;
661  unsigned int bandUsageIndex = 0;
662  unsigned int fileUsageIndex = 0;
663  unsigned int noBandsUsed = 0;
664  unsigned int localComparisonBandLim = 0;
665  bool toBeIgnored = false;
666  unsigned int maxShellsToUse = 0;
667  unsigned int minShellsToUse = 0;
668  bool firstDoneAlready = false;
669 
670  //======================================== Initialise internal values
671  std::vector<double> trSigmaWeights = std::vector<double> ( this->all->size() + 1 );
672  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( trSigmaWeights.size() ); iter++ )
673  {
674  trSigmaWeights.at(iter) = 0.0;
675  }
676 
677  //======================================== For each comparison against the one
678  fileUsageIndex = 0;
679  for ( unsigned int strIt = 0; strIt < static_cast<unsigned int> ( this->all->size() ); strIt++ )
680  {
681  //==================================== If not needed, do not compute
682  if ( this->_enLevelsDoNotFollow.at(strIt) == 1 ) { continue; }
683 
684  //==================================== Figure required number of bands
685  localComparisonBandLim = std::min( this->one->_bandwidthLimit, this->all->at(strIt)->_bandwidthLimit );
686  for ( unsigned int bandIter = 0; bandIter < localComparisonBandLim ; bandIter++ )
687  {
688  //================================ Ignore odd l's as they are 0 anyway
689  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
690 
691  //================================ Ignore l's designated to be ignored
692  toBeIgnored = false;
693  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
694  if ( toBeIgnored ) { continue; }
695 
696  //================================ This band will be used!
697  noBandsUsed += 1;
698  }
699 
700  //==================================== Allocate memory
701  bandUsageIndex = 0;
702  this->_EMatrices.emplace_back ( std::vector< std::vector< std::vector< std::array<double,2> > > > ( noBandsUsed ) );
703  for ( unsigned int bandIter = 0; bandIter < localComparisonBandLim; bandIter++ )
704  {
705  //================================ Ignore odd l's as they are 0 anyway
706  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
707 
708  //================================ Ignore l's designated to be ignored
709  toBeIgnored = false;
710  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
711  if ( toBeIgnored ) { continue; }
712 
713  //================================ Allocate vector of order1
714  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex) = std::vector< std::vector< std::array<double,2> > > ( static_cast<unsigned int> ( ( bandIter * 2 ) + 1 ) );
715 
716  //================================ Allocate vector of order2
717  for ( unsigned int lIter = 0; lIter < static_cast<unsigned int> ( ( bandIter * 2 ) + 1 ); lIter++ )
718  {
719  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(lIter) = std::vector< std::array<double,2> > ( static_cast<unsigned int> ( ( bandIter * 2 ) + 1 ) );
720  }
721 
722  bandUsageIndex += 1;
723  }
724 
725  //==================================== Set up
726  maxShellsToUse = std::max ( this->one->_noShellsWithData, this->all->at(strIt)->_noShellsWithData );
727  minShellsToUse = std::min ( this->one->_noShellsWithData, this->all->at(strIt)->_noShellsWithData );
728 
729  //==================================== Compute weights and abscissas for GL Integration
730  glIntegrationOrder = std::max ( glIntegOrderVec->at(0), glIntegOrderVec->at(strIt+1) );
731  glAbscissas = std::vector<double> ( glIntegrationOrder );
732  glWeights = std::vector<double> ( glIntegrationOrder );
733  ProSHADE_internal_legendre::getLegendreAbscAndWeights ( glIntegrationOrder, &glAbscissas, &glWeights, settings );
734 
735  //==================================== For each band (l)
736  bandUsageIndex = 0;
737  for ( unsigned int bandIter = 0; bandIter < localComparisonBandLim; bandIter++ )
738  {
739  //================================ Ignore odd l's as they are 0 anyway
740  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
741 
742  //================================ Ignore l's designated to be ignored
743  toBeIgnored = false;
744  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
745  if ( toBeIgnored ) { continue; }
746 
747  //================================ For each order (m)
748  for ( unsigned int orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
749  {
750  //============================ Clear values for weights
751  radiiVals1.clear ( );
752  radiiVals2.clear ( );
753  radiiVals1 = std::vector<double> ( this->one->_noShellsWithData );
754  radiiVals2 = std::vector<double> ( this->all->at(strIt)->_noShellsWithData );
755 
756  //============================ Get index for m
757  index11 = seanindex ( orderIter-bandIter, bandIter, this->one->_bandwidthLimit );
758  index12 = seanindex ( orderIter-bandIter, bandIter, this->all->at(strIt)->_bandwidthLimit );
759 
760  //============================ For each radius, deal with weights
761  for ( unsigned int radiusIter = 0; radiusIter < maxShellsToUse; radiusIter++ )
762  {
763  if ( !firstDoneAlready )
764  {
765  if ( radiusIter < this->one->_noShellsWithData )
766  {
767  //================ Multiply with itself - i.e. square
768  radiiVals1.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjugReal ( &this->one->_realSHCoeffs[radiusIter][index11],
769  &this->one->_imagSHCoeffs[radiusIter][index11],
770  &this->one->_realSHCoeffs[radiusIter][index11],
771  &this->one->_imagSHCoeffs[radiusIter][index11] );
772 
773  //================ Apply the integral r^2
774  radiiVals1.at(radiusIter) *= pow ( ( static_cast<double> ( radiusIter + 1 ) * shellSpacing ), 2.0 );
775  }
776  }
777 
778  if ( radiusIter < this->all->at(strIt)->_noShellsWithData )
779  {
780  //==================== Multiply with itself - i.e. square
781  radiiVals2.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjugReal ( &this->all->at(strIt)->_realSHCoeffs[radiusIter][index12],
782  &this->all->at(strIt)->_imagSHCoeffs[radiusIter][index12],
783  &this->all->at(strIt)->_realSHCoeffs[radiusIter][index12],
784  &this->all->at(strIt)->_imagSHCoeffs[radiusIter][index12] );
785 
786  //==================== Apply the integral r^2
787  radiiVals2.at(radiusIter) *= pow ( ( static_cast<double> ( radiusIter + 1 ) * shellSpacing ), 2.0 );
788  }
789  }
790 
791  //============================ Integrate weights
792  if ( !firstDoneAlready ) { trSigmaWeights.at(0) += this->glIntRR ( &radiiVals1, shellSpacing, &glAbscissas, &glWeights ); }
793  trSigmaWeights.at(strIt+1) += this->glIntRR ( &radiiVals2, shellSpacing, &glAbscissas, &glWeights );
794 
795  //============================ For each combination of m and m' for E matrices
796  for ( unsigned int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
797  {
798  //======================== Set index for order 2 (i.e. m' in E_{l,m,m'})
799  index21 = seanindex ( orderIter-bandIter, bandIter, this->one->_bandwidthLimit );
800  index22 = seanindex ( order2Iter-bandIter, bandIter, this->all->at(strIt)->_bandwidthLimit );
801 
802  //======================== Clear the vector to integrate over
803  radiiVals.clear ( );
804  radiiVals = std::vector< std::array<double,2> > ( minShellsToUse );
805 
806  //======================== Find the c*c values for different radii
807  for ( unsigned int radiusIter = 0; radiusIter < minShellsToUse; radiusIter++ )
808  {
809  //==================== Multiply coeffs
810  radiiVals.at(radiusIter) = ProSHADE_internal_misc::complexMultiplicationConjug ( &this->one->_realSHCoeffs[radiusIter][index21],
811  &this->one->_imagSHCoeffs[radiusIter][index21],
812  &this->all->at(strIt)->_realSHCoeffs[radiusIter][index22],
813  &this->all->at(strIt)->_imagSHCoeffs[radiusIter][index22] );
814 
815  //==================== Apply r^2 integral weight
816  radiiVals.at(radiusIter)[0] *= pow ( ( static_cast<double> ( radiusIter + 1 ) * shellSpacing ), 2.0 );
817  radiiVals.at(radiusIter)[1] *= pow ( ( static_cast<double> ( radiusIter + 1 ) * shellSpacing ), 2.0 );
818  }
819 
820  //======================== Integrate over all radii using 20-point Gauss-Legendre
821  arrVal1 = this->glIntCR ( &radiiVals, shellSpacing, &glAbscissas, &glWeights );
822 
823  //======================== Save the result into E matrices
824  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(orderIter).at(order2Iter)[0] = arrVal1[0];
825  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(orderIter).at(order2Iter)[1] = arrVal1[1];
826  }
827  }
828  bandUsageIndex += 1;
829  }
830 
831  //==================================== For each band (l), normalise using the correlation-like formula
832  bandUsageIndex = 0;
833  correlNormFactor = sqrt ( trSigmaWeights.at(0) * trSigmaWeights.at(strIt+1) );
834 
835  for ( unsigned int bandIter = 0; bandIter < localComparisonBandLim; bandIter++ )
836  {
837  //================================ Ignore odd l's as they are 0 anyway
838  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
839 
840  //================================ Ignore l's designated to be ignored
841  toBeIgnored = false;
842  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
843  if ( toBeIgnored ) { continue; }
844 
845  //================================ For each combination of m and m' for E matrices
846  for ( unsigned int orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
847  {
848  for ( unsigned int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
849  {
850  //======================== Normalisation using the correlation formula
851  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(orderIter).at(order2Iter)[0] /= correlNormFactor;
852  this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(orderIter).at(order2Iter)[1] /= correlNormFactor;
853  }
854  }
855  bandUsageIndex += 1;
856  }
857 
858  firstDoneAlready = true;
859  fileUsageIndex += 1;
860  }
861 
862  //======================================== Done
863  this->_trSigmaPreComputed = true;
864 
865 }
866 
879 std::vector<double> ProSHADE_internal::ProSHADE_compareOneAgainstAll::getTrSigmaDistance ( int verbose, ProSHADE::ProSHADE_settings* settings )
880 {
881  //======================================== Sanity check
882  if ( !this->_trSigmaPreComputed )
883  {
884  std::cerr << "!!! ProSHADE ERROR !!! Error in file compatison !!! Attempted to get TrSigma distance without pre-computing the required data (E matrices). Call the precomputeTrSigmaDescriptor function BEFORE calling the getTrSigmaDistance function." << std::endl;
885 
886  if ( settings->htmlReport )
887  {
888  std::stringstream hlpSS;
889  hlpSS << "<font color=\"red\">" << "Attempted to get TrSigma distance without pre-computing the required data (E matrices). This looks like an internal bug, please report this case." << "</font>";
890  rvapi_set_text ( hlpSS.str().c_str(),
891  "ProgressSection",
892  settings->htmlReportLineProgress,
893  1,
894  1,
895  1 );
896  settings->htmlReportLineProgress += 1;
897  rvapi_flush ( );
898  }
899 
900  exit ( -1 );
901  }
902 
903  //======================================== If already done, do not repeat :-)
904  if ( this->_trSigmaComputed )
905  {
906  return ( this->_distancesTraceSigma );
907  }
908  if ( verbose > 3 )
909  {
910  std::cout << ">>>>>>>> Sanity checks passed." << std::endl;
911  }
912 
913  //======================================== Initialise local variables
914  this->_distancesTraceSigma = std::vector<double> ( this->all->size() );
915  std::vector<double> hlpVec;
916  unsigned int bandUsageIndex = 0;
917  bool toBeIgnored = false;
918  unsigned int localComparisonBandLim = 0;
919  unsigned int fileUsageIndex = 0;
920 
921  //======================================== For each comparison against the one
922  for ( unsigned int strIt = 0; strIt < static_cast<unsigned int> ( this->all->size() ); strIt++ )
923  {
924  //===================================== If not needed, do not compute
925  if ( this->_enLevelsDoNotFollow.at(strIt) == 1 ) { this->_distancesTraceSigma.at(strIt) = -999.9; continue; }
926 
927  if ( verbose > 3 )
928  {
929  std::cout << ">>>>>>>> Now computing distance between structure 0 and structure " << strIt+1 << " ." << std::endl;
930  }
931 
932  //==================================== Set up
933  bandUsageIndex = 0;
934  localComparisonBandLim = std::min( this->one->_bandwidthLimit, this->all->at(strIt)->_bandwidthLimit );
935  this->_distancesTraceSigma.at(strIt) = 0.0;
936 
937  //==================================== For each l (band)
938  for ( unsigned int lIter = 0; lIter < localComparisonBandLim; lIter++ )
939  {
940  //================================ Ignore odd l's as they are 0 anyway
941  if ( !this->_keepOrRemove ) { if ( ( lIter % 2 ) != 0 ) { continue; } }
942 
943  //================================ Ignore l's designated to be ignored
944  toBeIgnored = false;
945  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( lIter == static_cast<unsigned int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
946  if ( toBeIgnored ) { continue; }
947 
948  //================================ Find the complex matrix SVD singular values
949  hlpVec = this->getSingularValues ( fileUsageIndex, bandUsageIndex, static_cast<unsigned int> ( ( lIter * 2 ) + 1 ), settings );
950 
951  //================================ Now sum the trace
952  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( hlpVec.size() ); iter++ )
953  {
954  this->_distancesTraceSigma.at(strIt) += hlpVec.at(iter);
955  }
956 
957  bandUsageIndex += 1;
958  }
959 
960  fileUsageIndex += 1;
961 
962  if ( verbose > 2 )
963  {
964  std::cout << ">>>>> Distance between structure 0 and structure " << strIt+1 << " computed." << std::endl;
965  }
966  }
967 
968  //======================================== Done
969  this->_trSigmaComputed = true;
970 
971  //======================================== Return
972  return ( this->_distancesTraceSigma );
973 
974 }
975 
988 double ProSHADE_internal::ProSHADE_comparePairwise::getRotCoeffDistance ( ProSHADE::ProSHADE_settings* settings )
989 {
990  //======================================== Sanity check
991  if ( !this->_wignerMatricesComputed )
992  {
993  std::cerr << "!!! ProSHADE ERROR !!! Error in file compatison !!! Attempted to get Full Rotation distance without pre-computing the Wigner D matrices. Call the generateWignerMatrices function BEFORE calling the getRotCoeffDistance function." << std::endl;
994 
995  if ( settings->htmlReport )
996  {
997  std::stringstream hlpSS;
998  hlpSS << "<font color=\"red\">" << "Attempted to get Full Rotation distance without pre-computing the Wigner D matrices. This looks like an internal bug, please report this case." << "</font>";
999  rvapi_set_text ( hlpSS.str().c_str(),
1000  "ProgressSection",
1001  settings->htmlReportLineProgress,
1002  1,
1003  1,
1004  1 );
1005  settings->htmlReportLineProgress += 1;
1006  rvapi_flush ( );
1007  }
1008 
1009  exit ( -1 );
1010  }
1011 
1012  //======================================== Initialise local variables
1013  bool toBeIgnored = false;
1014  double totSum = 0.0;
1015  unsigned int bandUsageIndex = 0;
1016 
1017  //======================================== For each band
1018  for ( int bandIter = 0; bandIter < static_cast<int> ( this->_bandwidthLimit ); bandIter++ )
1019  {
1020  //==================================== Ignore odd l's as they are 0 anyway
1021  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
1022 
1023  //==================================== Ignore l's designated to be ignored
1024  toBeIgnored = false;
1025  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
1026  if ( toBeIgnored ) { continue; }
1027 
1028  //==================================== For each order1
1029  for ( int order1 = 0; order1 < static_cast<int> ( ( bandIter * 2 ) + 1 ); order1++ )
1030  {
1031  //================================ For each order2
1032  for ( int order2 = 0; order2 < static_cast<int> ( ( bandIter * 2 ) + 1 ); order2++ )
1033  {
1034  //============================ Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
1035  totSum += ProSHADE_internal_misc::complexMultiplicationReal ( &this->_wignerMatrices.at(bandIter).at(order1).at(order2)[0],
1036  &this->_wignerMatrices.at(bandIter).at(order1).at(order2)[1],
1037  &this->_trSigmaEMatrix[bandUsageIndex][order2][order1][0],
1038  &this->_trSigmaEMatrix[bandUsageIndex][order2][order1][1] );
1039  }
1040  }
1041 
1042  //==================================== Prepare for next loop
1043  bandUsageIndex += 1;
1044  }
1045 
1046  //======================================== Done
1047  return ( totSum );
1048 }
1049 
1063 std::vector<double> ProSHADE_internal::ProSHADE_compareOneAgainstAll::getRotCoeffDistance ( int verbose, ProSHADE::ProSHADE_settings* settings )
1064 {
1065  //======================================== Sanity check
1066  if ( !this->_wignerMatricesComputed )
1067  {
1068  std::cerr << "!!! ProSHADE ERROR !!! Error in file compatison !!! Attempted to get Full Rotation distance without pre-computing the Wigner D matrices. Call the generateWignerMatrices function BEFORE calling the getRotCoeffDistance function." << std::endl;
1069 
1070  if ( settings->htmlReport )
1071  {
1072  std::stringstream hlpSS;
1073  hlpSS << "<font color=\"red\">" << "Attempted to get Full Rotation distance without pre-computing the Wigner D matrices. This looks like an internal bug, please report this case." << "</font>";
1074  rvapi_set_text ( hlpSS.str().c_str(),
1075  "ProgressSection",
1076  settings->htmlReportLineProgress,
1077  1,
1078  1,
1079  1 );
1080  settings->htmlReportLineProgress += 1;
1081  rvapi_flush ( );
1082  }
1083 
1084  exit ( -1 );
1085  }
1086 
1087  //======================================== Save some time on re-runs
1088  if ( this->_fullDistComputed )
1089  {
1090  return ( this->_distancesFullRotation );
1091  }
1092  if ( verbose > 3 )
1093  {
1094  std::cout << ">>>>>>>> Sanity checks passed." << std::endl;
1095  }
1096 
1097  //======================================== Initialise local variables
1098  bool toBeIgnored = false;
1099  double totSum = 0.0;
1100  unsigned int bandUsageIndex = 0;
1101  unsigned int fileUsageIndex = 0;
1102  int localComparisonBandLim = 0;
1103 
1104  for ( unsigned int strIt = 0; strIt < static_cast<unsigned int> ( this->all->size() ); strIt++ )
1105  {
1106  //==================================== If not needed, do not compute
1107  if ( this->_enLevelsDoNotFollow.at(strIt) == 1 ) { this->_distancesFullRotation.emplace_back( -999.9 ); continue; }
1108  if ( this->_trSigmaDoNotFollow.at(strIt) == 1 ) { this->_distancesFullRotation.emplace_back( -999.9 ); fileUsageIndex += 1; continue; }
1109 
1110  if ( verbose > 3 )
1111  {
1112  std::cout << ">>>>>>>> Now computing distance between structure 0 and structure " << strIt+1 << " ." << std::endl;
1113  }
1114 
1115  //==================================== Find the appropriate bandwidth
1116  localComparisonBandLim = std::min ( this->one->_bandwidthLimit, this->all->at(strIt)->_bandwidthLimit );
1117 
1118  //==================================== Initialise local variables
1119  toBeIgnored = false;
1120  totSum = 0.0;
1121  bandUsageIndex = 0;
1122 
1123  //==================================== For each band
1124  for ( int bandIter = 0; bandIter < static_cast<int> ( localComparisonBandLim ); bandIter++ )
1125  {
1126  //================================ Ignore odd l's as they are 0 anyway
1127  if ( !this->_keepOrRemove ) { if ( ( bandIter % 2 ) != 0 ) { continue; } }
1128 
1129  //================================ Ignore l's designated to be ignored
1130  toBeIgnored = false;
1131  for ( unsigned int igIt = 0; igIt < static_cast<unsigned int> ( this->_lsToIgnore.size() ); igIt++ ) { if ( bandIter == static_cast<int> ( this->_lsToIgnore.at(igIt) ) ) { toBeIgnored = true; } }
1132  if ( toBeIgnored ) { continue; }
1133 
1134  //================================ For each order1
1135  for ( int order1 = 0; order1 < static_cast<int> ( ( bandIter * 2 ) + 1 ); order1++ )
1136  {
1137  //============================ For each order2
1138  for ( int order2 = 0; order2 < static_cast<int> ( ( bandIter * 2 ) + 1 ); order2++ )
1139  {
1140  //======================== Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
1141  totSum += ProSHADE_internal_misc::complexMultiplicationReal ( &this->_wignerMatrices.at(strIt).at(bandIter).at(order1).at(order2)[0],
1142  &this->_wignerMatrices.at(strIt).at(bandIter).at(order1).at(order2)[1],
1143  &this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(order2).at(order1)[0],
1144  &this->_EMatrices.at(fileUsageIndex).at(bandUsageIndex).at(order2).at(order1)[1] );
1145  }
1146  }
1147 
1148  //================================ Prepare for next loop
1149  bandUsageIndex += 1;
1150  }
1151 
1152  //==================================== Save
1153  this->_distancesFullRotation.emplace_back ( totSum );
1154  fileUsageIndex += 1;
1155 
1156  if ( verbose > 2 )
1157  {
1158  std::cout << ">>>>> Distance between structure 0 and structure " << strIt+1 << " computed." << std::endl;
1159  }
1160  }
1161 
1162  //======================================== Done
1163  this->_fullDistComputed = true;
1164 
1165  //======================================== Return
1166  return ( this->_distancesFullRotation );
1167 }
ProSHADE_internal::ProSHADE_compareOneAgainstAll::getRotCoeffDistance
std::vector< double > getRotCoeffDistance(int verbose, ProSHADE::ProSHADE_settings *settings)
This function computes the full rotation function descriptor distances.
Definition: ProSHADE_descriptors.cpp:1063
ProSHADE::ProSHADE_settings::htmlReport
bool htmlReport
Should HTML report for the run be created?
Definition: ProSHADE.h:186
ProSHADE_internal::ProSHADE_compareOneAgainstAll::getTrSigmaDistance
std::vector< double > getTrSigmaDistance(int verbose, ProSHADE::ProSHADE_settings *settings)
This function computes the trace sigma descriptor distances.
Definition: ProSHADE_descriptors.cpp:879
ProSHADE.h
The main header file containing all declarations the user of the library needs.
ProSHADE_internal::ProSHADE_comparePairwise::getRotCoeffDistance
double getRotCoeffDistance(ProSHADE::ProSHADE_settings *settings)
This function computes the full rotation function descriptor distances.
Definition: ProSHADE_descriptors.cpp:988
ProSHADE_internal::ProSHADE_comparePairwise::precomputeTrSigmaDescriptor
void precomputeTrSigmaDescriptor()
This function computes the E matrices required for the trace sigma descriptor, the rotation function ...
Definition: ProSHADE_descriptors.cpp:454
ProSHADE_legendre.h
This file contains all the functions related to computing the Gauss-Legendre integration variables...
ProSHADE_internal_misc::complexMultiplication
std::array< double, 2 > complexMultiplication(double *r1, double *i1, double *r2, double *i2)
Function to multiply two complex numbers.
Definition: ProSHADE_descriptors.cpp:112
ProSHADE_internal::ProSHADE_compareOneAgainstAll::precomputeTrSigmaDescriptor
void precomputeTrSigmaDescriptor(double shellSpacing, std::vector< unsigned int > *glIntegOrderVec, ProSHADE::ProSHADE_settings *settings)
This function computes the E matrices required for the trace sigma descriptor, the rotation function ...
Definition: ProSHADE_descriptors.cpp:646
ProSHADE_internal.h
The main header file containing all declarations for the innter workings of the library.
ProSHADE::ProSHADE_settings::htmlReportLineProgress
int htmlReportLineProgress
Iterator for current HTML line in the progress bar.
Definition: ProSHADE.h:188
ProSHADE_internal_misc::complexMultiplicationConjugReal
double complexMultiplicationConjugReal(double *r1, double *i1, double *r2, double *i2)
Function to multiply two complex numbers where the second number is conjugated before the multiplicat...
Definition: ProSHADE_descriptors.cpp:199
ProSHADE_internal::ProSHADE_compareOneAgainstAll::getEnergyLevelsDistance
std::vector< double > getEnergyLevelsDistance(int verbose, ProSHADE::ProSHADE_settings *settings)
This function computes the energy level descriptor value from the first structure to all remaining st...
Definition: ProSHADE_descriptors.cpp:310
ProSHADE_internal::ProSHADE_data::precomputeRotInvDescriptor
void precomputeRotInvDescriptor(ProSHADE::ProSHADE_settings *settings)
This function computes the RRP matrices, which are required for the computation of the energy levels ...
Definition: ProSHADE_descriptors.cpp:219
ProSHADE::ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of multitude of par...
Definition: ProSHADE.h:74
ProSHADE_internal_misc::complexMultiplicationConjug
std::array< double, 2 > complexMultiplicationConjug(double *r1, double *i1, double *r2, double *i2)
Function to multiply two complex numbers where the second number is conjugated before the multiplicat...
Definition: ProSHADE_descriptors.cpp:143
ProSHADE_misc.h
This file contains the ProSHADE_internal_misc namespace and its miscellaneous functions.
ProSHADE_internal_misc::pearsonCorrCoeff
double pearsonCorrCoeff(double *valSet1, double *valSet2, unsigned int length)
Function for computing the Pearson correlation coefficient.
Definition: ProSHADE_descriptors.cpp:62
ProSHADE_internal_misc::complexMultiplicationReal
double complexMultiplicationReal(double *r1, double *i1, double *r2, double *i2)
Function to multiply two complex numbers, returning only the real part of the result.
Definition: ProSHADE_descriptors.cpp:174