ProSHADE  0.6.5 (NOV 2018)
Protein Shape Descriptors and Symmetry Detection
ProSHADE_lapack.cpp
Go to the documentation of this file.
1 
20 //============================================ Lapack-E
21 #include <complex>
22 #define lapack_complex_double std::complex<double>
23 #define lapack_complex_float std::complex<float>
24 #define LAPACK_COMPLEX_CUSTOM
25 #include <lapacke.h>
26 
27 //============================================ RVAPI
28 #include <rvapi_interface.h>
29 
30 //============================================ ProSHADE
31 #include "ProSHADE.h"
32 #include "ProSHADE_internal.h"
33 
47 std::vector<double> ProSHADE_internal::ProSHADE_comparePairwise::getSingularValues ( unsigned int band,
48  unsigned int dim,
49  ProSHADE::ProSHADE_settings* settings )
50 {
51  //======================================== Initialise local variables
52  double *singularValues = new double[dim];
53  __extension__ std::complex<double> rotMatU[dim][dim];
54  __extension__ std::complex<double> rotMatV[dim][dim];
55  std::vector<double> ret;
56  int returnValue = 0;
57 
58  //======================================== Load SH coeffs into array in column-major order
59  std::complex<double> *matrixToDecompose = new std::complex<double>[dim*dim];
60  for ( unsigned int rowIt = 0; rowIt < dim; rowIt++ )
61  {
62  for ( unsigned int colIt = 0; colIt < dim; colIt++ )
63  {
64  matrixToDecompose[(colIt*dim)+rowIt] = std::complex<double> ( this->_trSigmaEMatrix[band][rowIt][colIt][0],
65  this->_trSigmaEMatrix[band][rowIt][colIt][1] );
66  }
67  }
68 
69  //======================================== Call the LAPACK ZGESSD function
70  returnValue = LAPACKE_zgesdd ( LAPACK_COL_MAJOR, 'A', dim, dim, matrixToDecompose, dim, singularValues, *rotMatU, dim, *rotMatV, dim );
71 
72  //======================================== Check for errors
73  if ( returnValue != 0 )
74  {
75  std::cerr << "!!! ProSHADE ERROR !!! Complex Matrix SVD algorithm did not converge, cannot proceed!" << std::endl;
76 
77  if ( settings->htmlReport )
78  {
79  std::stringstream hlpSS;
80  hlpSS << "<font color=\"red\">" << "The LAPACK Complex Matrix SVD algorithm did not converge, cannot proceed with computation of Trace Sigma distances. This problem could be resolved by changing the resolution or not using the Trace Sigma distances; these are not good solutions, but they are all I have got right now." << "</font>";
81  rvapi_set_text ( hlpSS.str().c_str(),
82  "ProgressSection",
83  settings->htmlReportLineProgress,
84  1,
85  1,
86  1 );
87  settings->htmlReportLineProgress += 1;
88  rvapi_flush ( );
89  }
90 
91  exit ( -1 );
92  }
93 
94  //======================================== Save results in output format
95  for ( unsigned int iter = 0; iter < dim; iter++ )
96  {
97  ret.emplace_back ( singularValues[iter] );
98  }
99 
100  //======================================== Free memory
101  delete[] singularValues;
102  delete[] matrixToDecompose;
103 
104  //======================================== Done
105  return ( ret );
106 
107 }
108 
122 std::vector<double> ProSHADE_internal::ProSHADE_compareOneAgainstAll::getSingularValues ( unsigned int strNo, unsigned int band, unsigned int dim, ProSHADE::ProSHADE_settings* settings )
123 {
124  //======================================== Initialise local variables
125  double *singularValues = new double[dim];
126  __extension__ std::complex<double> rotMatU[dim][dim];
127  __extension__ std::complex<double> rotMatV[dim][dim];
128  std::vector<double> ret;
129  int returnValue = 0;
130 
131  //======================================== Load SH coeffs into array in column-major order
132  std::complex<double> *matrixToDecompose = new std::complex<double>[dim*dim];
133  for ( unsigned int rowIt = 0; rowIt < dim; rowIt++ )
134  {
135  for ( unsigned int colIt = 0; colIt < dim; colIt++ )
136  {
137  matrixToDecompose[(colIt*dim)+rowIt] = std::complex<double> ( this->_EMatrices.at(strNo).at(band).at(rowIt).at(colIt)[0],
138  this->_EMatrices.at(strNo).at(band).at(rowIt).at(colIt)[1] );
139  }
140  }
141 
142  //======================================== Call the LAPACK ZGESSD routine
143  returnValue = LAPACKE_zgesdd ( LAPACK_COL_MAJOR, 'A', dim, dim, matrixToDecompose, dim, singularValues, *rotMatU, dim, *rotMatV, dim );
144 
145  //======================================== Check for errors
146  if ( returnValue != 0 )
147  {
148  std::cerr << "!!! ProSHADE ERROR !!! Complex Matrix SVD algorithm did not converge, cannot proceed!" << std::endl;
149 
150  if ( settings->htmlReport )
151  {
152  std::stringstream hlpSS;
153  hlpSS << "<font color=\"red\">" << "The LAPACK Complex Matrix SVD algorithm did not converge, cannot proceed with computation of Trace Sigma distances. This problem could be resolved by changing the resolution or not using the Trace Sigma distances; these are not good solutions, but they are all I have got right now." << "</font>";
154  rvapi_set_text ( hlpSS.str().c_str(),
155  "ProgressSection",
156  settings->htmlReportLineProgress,
157  1,
158  1,
159  1 );
160  settings->htmlReportLineProgress += 1;
161  rvapi_flush ( );
162  }
163 
164  exit ( -1 );
165  }
166 
167  //======================================== Save results in output format
168  for ( unsigned int iter = 0; iter < dim; iter++ )
169  {
170  ret.emplace_back ( singularValues[iter] );
171  }
172 
173  //======================================== Free memory
174  delete[] singularValues;
175  delete[] matrixToDecompose;
176 
177  //======================================== Done
178  return ( ret );
179 
180 }
181 
194 std::vector< std::vector< std::vector<double> > > ProSHADE_internal::ProSHADE_comparePairwise::getSingularValuesUandVMatrix ( std::vector< std::vector<double> > mat,
195  unsigned int dim,
196  ProSHADE::ProSHADE_settings* settings )
197 {
198  //======================================== Initialise local variables
199  double *singularValues = new double[dim];
200  __extension__ std::complex<double> rotMatU[dim][dim];
201  __extension__ std::complex<double> rotMatV[dim][dim];
202  std::vector< std::vector< std::vector<double> > > ret;
203  std::vector< std::vector<double> > retHlp;
204  std::vector<double> hlpVec;
205  int returnValue = 0;
206 
207  //======================================== Load SH coeffs into array in column-major order
208  std::complex<double> *matrixToDecompose = new std::complex<double>[dim*dim];
209  for ( unsigned int rowIt = 0; rowIt < dim; rowIt++ )
210  {
211  for ( unsigned int colIt = 0; colIt < dim; colIt++ )
212  {
213  matrixToDecompose[(colIt*dim)+rowIt] = std::complex<double> ( mat.at(rowIt).at(colIt),
214  0.0 );
215  }
216  }
217 
218  //======================================== Call the LAPACK ZGESSD function
219  returnValue = LAPACKE_zgesdd ( LAPACK_COL_MAJOR, 'A', dim, dim, matrixToDecompose, dim, singularValues, *rotMatU, dim, *rotMatV, dim );
220 
221  //======================================== Check for errors
222  if ( returnValue != 0 )
223  {
224  std::cerr << "!!! ProSHADE ERROR !!! Complex Matrix SVD algorithm did not converge, cannot proceed!" << std::endl;
225 
226  if ( settings->htmlReport )
227  {
228  std::stringstream hlpSS;
229  hlpSS << "<font color=\"red\">" << "The LAPACK Complex Matrix SVD algorithm did not converge, cannot proceed with computation of Trace Sigma distances. This problem could be resolved by changing the resolution or not using the Trace Sigma distances; these are not good solutions, but they are all I have got right now." << "</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  //======================================== Save results in output format
244  for ( unsigned int iter = 0; iter < dim; iter++ ) { hlpVec.emplace_back ( 0.0 ); }
245  for ( unsigned int iter = 0; iter < dim; iter++ ) { retHlp.emplace_back ( hlpVec ); }
246  ret.emplace_back ( retHlp ); ret.emplace_back ( retHlp );
247 
248  //======================================== Save U
249  for ( unsigned int rowIt = 0; rowIt < dim; rowIt++ )
250  {
251  for ( unsigned int colIt = 0; colIt < dim; colIt++ )
252  {
253  ret.at(0).at(rowIt).at(colIt) = rotMatU[rowIt][colIt].real();
254  }
255  }
256 
257  //======================================== Save V
258  for ( unsigned int rowIt = 0; rowIt < dim; rowIt++ )
259  {
260  for ( unsigned int colIt = 0; colIt < dim; colIt++ )
261  {
262  ret.at(1).at(rowIt).at(colIt) = rotMatV[rowIt][colIt].real();
263  }
264  }
265 
266  //======================================== Free memory
267  delete[] singularValues;
268  delete[] matrixToDecompose;
269 
270  //======================================== Done
271  return ( ret );
272 
273 }
ProSHADE::ProSHADE_settings::htmlReport
bool htmlReport
Should HTML report for the run be created?
Definition: ProSHADE.h:186
ProSHADE.h
The main header file containing all declarations the user of the library needs.
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::ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of multitude of par...
Definition: ProSHADE.h:74