ʻOhana
Population structure, admixture history, and selection using learning methods.
jade.matrix.hpp
1 /* -------------------------------------------------------------------------
2  Ohana
3  Copyright (c) 2015-2020 Jade Cheng (\___/)
4  Jade Cheng <info@jade-cheng.com> (='.'=)
5  ------------------------------------------------------------------------- */
6 
7 #ifndef JADE_MATRIX_HPP__
8 #define JADE_MATRIX_HPP__
9 
10 #include "jade.blas.hpp"
11 #include "jade.error.hpp"
12 #include "jade.lapack.hpp"
13 
14 namespace jade
15 {
16  ///
17  /// A template class for a matrix.
18  ///
19  template <typename TValue>
21  {
22  public:
23  /// The value type.
24  typedef TValue value_type;
25 
26  /// The BLAS type.
28 
29  /// The LAPACK type.
31 
32  /// The initializer list type.
33  typedef std::initializer_list<
34  std::initializer_list<value_type>
36 
37  ///
38  /// Initializes a new instance of the class with no size.
39  ///
40  inline basic_matrix()
41  : basic_matrix(0, 0)
42  {
43  }
44 
45  ///
46  /// Initializes a new instance of the class based on the specified
47  /// width and height.
48  ///
50  const size_t cy, ///< The height.
51  const size_t cx) ///< The width.
52  : _cy (cy)
53  , _cx (cx)
54  , _m ()
55  {
56  assert((cx == 0 && cy == 0) || (cx != 0 && cy != 0));
57  _m.resize(cx * cy);
58  }
59 
60  ///
61  /// Initializes a new instance of the class based on values from the
62  /// specified file.
63  ///
64  inline explicit basic_matrix(
65  const char * const path) ///< The path to the file.
66  : basic_matrix()
67  {
68  read(path);
69  }
70 
71  ///
72  /// Initializes a new instance of the class based on values from the
73  /// specified file.
74  ///
75  inline explicit basic_matrix(
76  const std::string & path) ///< The path to the file.
77  : basic_matrix(path.c_str())
78  {
79  }
80 
81  ///
82  /// Initializes a new instance of the class based on values from the
83  /// specified input stream.
84  ///
85  inline explicit basic_matrix(
86  std::istream & in) ///< The input stream.
87  : basic_matrix()
88  {
89  read(in);
90  }
91 
92  ///
93  /// Initializes a new instance of the class based on the specified
94  /// values.
95  ///
96  /// \param values The two-dimensional values.
97  ///
98  explicit basic_matrix(const initializer_list_type & values)
99  : _cy (values.size())
100  , _cx (_cy == 0 ? 0 : values.begin()->size())
101  , _m ()
102  {
103  assert(std::all_of(values.begin(), values.end(), [this](
104  const std::initializer_list<value_type> & row_values) -> bool {
105  return row_values.size() == _cx; }));
106 
107  _m.reserve(_cy * _cx);
108  for (const auto & row : values)
109  for (const auto n : row)
110  _m.push_back(n);
111  }
112 
113  ///
114  /// \return True if the specified predicate is true for all values.
115  ///
116  template <typename TPredicate>
117  inline bool all_of(
118  const TPredicate predicate) ///< The predicate function.
119  const
120  {
121  return std::all_of(_m.begin(), _m.end(), predicate);
122  }
123 
124  ///
125  /// \return True if the specified predicate is true for any value.
126  ///
127  template <typename TPredicate>
128  inline bool any_of(
129  const TPredicate predicate) ///< The predicate function.
130  const
131  {
132  return std::any_of(_m.begin(), _m.end(), predicate);
133  }
134 
135  ///
136  /// Clamps all values in the matrix to the specified range.
137  ///
138  void clamp(
139  const value_type min, ///< The minimum value.
140  const value_type max) ///< The maximum value.
141  {
142  assert(min <= max);
143 
144  auto dst_ptr = get_data();
145  const auto dst_end = dst_ptr + get_length();
146 
147  while (dst_ptr != dst_end)
148  {
149  *dst_ptr = std::min(std::max(min, *dst_ptr), max);
150  dst_ptr++;
151  }
152  }
153 
154  ///
155  /// Clamps all values in a column to the specified range.
156  ///
158  const size_t column, ///< The column to clamp.
159  const value_type min, ///< The minimum value.
160  const value_type max) ///< The maximum value.
161  {
162  assert(column < get_width());
163  assert(min <= max);
164 
165  auto dst_ptr = get_data(0, column);
166  const auto dst_end = dst_ptr + get_length();
167 
168  while (dst_ptr != dst_end)
169  {
170  *dst_ptr = std::min(std::max(min, *dst_ptr), max);
171  dst_ptr += get_width();
172  }
173  }
174 
175  ///
176  /// Clamps all values in a row to the specified range.
177  ///
178  void clamp_row(
179  const size_t row, ///< The row to clamp.
180  const value_type min, ///< The minimum value.
181  const value_type max) ///< The maximum value.
182  {
183  assert(row < get_height());
184  assert(min <= max);
185 
186  auto dst_ptr = get_data(row, 0);
187  const auto dst_end = dst_ptr + get_width();
188 
189  while (dst_ptr != dst_end)
190  {
191  *dst_ptr = std::min(std::max(min, *dst_ptr), max);
192  dst_ptr++;
193  }
194  }
195 
196  ///
197  /// \return True if any element in the matrix is infinite.
198  ///
199  bool contains_inf() const
200  {
201  auto src_ptr = get_data();
202  const auto src_end = src_ptr + get_length();
203 
204  while (src_ptr != src_end)
205  if (std::isinf(*src_ptr++))
206  return true;
207 
208  return false;
209  }
210 
211  ///
212  /// \return True if any element in the matrix is not a number.
213  ///
214  bool contains_nan() const
215  {
216  auto src_ptr = get_data();
217  const auto src_end = src_ptr + get_length();
218 
219  while (src_ptr != src_end)
220  if (std::isnan(*src_ptr++))
221  return true;
222 
223  return false;
224  }
225 
226  ///
227  /// \return A new column vector based on the values from the specified
228  /// column.
229  ///
231  const size_t column) ///< The column to copy.
232  const
233  {
234  basic_matrix out;
235  copy_column(column, out);
236  return out;
237  }
238 
239  ///
240  /// Copies a column into the specified vector.
241  ///
243  const size_t column, ///< The column to copy.
244  basic_matrix & out) ///< The output vector.
245  const
246  {
247  assert(column < get_width());
248 
249  out.resize(get_height(), 1);
250 
251  auto dst_ptr = out.get_data();
252  auto src_ptr = get_data(0, column);
253  const auto src_end = src_ptr + get_length();
254 
255  while (src_ptr != src_end)
256  {
257  *dst_ptr++ = *src_ptr;
258  src_ptr += get_width();
259  }
260  }
261 
262  ///
263  /// Copies the elements in the lower triangle of the square matrix to
264  /// the corresponding upper triangle elements.
265  ///
267  {
268  assert(is_square());
269 
270  const auto n = get_height();
271  if (n < 2)
272  return;
273 
274  //
275  // Prepare pointers to the data.
276  //
277  const auto dst_end = get_data() + get_length();
278  auto dst_ptr = get_data() + 1;
279  auto src_end = get_data() + get_length();
280  auto src_ptr = get_data() + n;
281 
282  //
283  // Loop until the entire triangle is copied.
284  //
285  while (dst_ptr != dst_end)
286  {
287  const auto dst_anchor = dst_ptr;
288  const auto src_anchor = src_ptr;
289 
290  //
291  // Advance down the row of the lower triangle and across the
292  // column of the upper triangle.
293  //
294  while (src_ptr != src_end)
295  {
296  *dst_ptr++ = *src_ptr;
297  src_ptr += n;
298  }
299 
300  //
301  // Adjust the pointers for the next row and column.
302  //
303  dst_ptr = dst_anchor + n + 1;
304  src_ptr = src_anchor + n + 1;
305  src_end++;
306  }
307  }
308 
309  ///
310  /// \return A new row vector based on the values from the specified row.
311  ///
313  const size_t row) ///< The row to copy.
314  const
315  {
316  basic_matrix out;
317  copy_row(row, out);
318  return out;
319  }
320 
321  ///
322  /// Copies a row into the specified matrix.
323  ///
324  void copy_row(
325  const size_t row, ///< The row to copy.
326  basic_matrix & out) ///< The output vector.
327  const
328  {
329  assert(row < get_height());
330 
331  out.resize(1, get_width());
332 
333  auto dst_ptr = out.get_data();
334  auto src_ptr = get_data(row, 0);
335  const auto src_end = src_ptr + get_width();
336 
337  while (src_ptr != src_end)
338  *dst_ptr++ = *src_ptr++;
339  }
340 
341  ///
342  /// Copies the elements in the upper triangle of the square matrix to
343  /// the corresponding lower triangle elements.
344  ///
346  {
347  assert(is_square());
348 
349  const auto n = get_height();
350  if (n < 2)
351  return;
352 
353  //
354  // Prepare pointers to the data.
355  //
356  const auto src_end = get_data() + get_length();
357  auto src_ptr = get_data() + 1;
358  auto dst_end = get_data() + get_length();
359  auto dst_ptr = get_data() + n;
360 
361  //
362  // Loop until the entire triangle is copied.
363  //
364  while (src_ptr != src_end)
365  {
366  const auto dst_anchor = dst_ptr;
367  const auto src_anchor = src_ptr;
368 
369  //
370  // Advance down the row of the lower triangle and across the
371  // column of the upper triangle.
372  //
373  while (dst_ptr != dst_end)
374  {
375  *dst_ptr = *src_ptr++;
376  dst_ptr += n;
377  }
378 
379  //
380  // Adjust the pointers for the next row and column.
381  //
382  dst_ptr = dst_anchor + n + 1;
383  src_ptr = src_anchor + n + 1;
384  dst_end++;
385  }
386  }
387 
388  ///
389  /// \return The transpose of the matrix.
390  ///
392  {
393  basic_matrix out;
394  create_transpose(out);
395  return out;
396  }
397 
398  ///
399  /// Stores the transpose of the matrix.
400  ///
402  basic_matrix & out) ///< The matrix receiving the transpose.
403  const
404  {
405  if (this == &out)
406  {
407  if (is_empty())
408  return;
409 
410  if (is_vector())
411  {
412  std::swap(out._cx, out._cy);
413  return;
414  }
415 
416  if (is_square())
417  {
418  for (size_t i = 1; i < _cy; i++)
419  for (size_t j = 0; j < i; j++)
420  std::swap(out(i, j), out(j, i));
421  return;
422  }
423 
424  basic_matrix tmp;
425  create_transpose(tmp);
426  out.swap(tmp);
427  return;
428  }
429 
430  if (is_empty() || is_vector() || is_square())
431  {
432  out = *this;
433  out.create_transpose(out);
434  return;
435  }
436 
437  out.resize(_cx, _cy);
438 
439  auto & src = *this;
440  for (size_t i = 0; i < _cy; i++)
441  for (size_t j = 0; j < _cx; j++)
442  out(j, i) = src(i, j);
443  }
444 
445  ///
446  /// Multiplies a left matrix by a right matrix and stores the result
447  /// into a destination matrix. If the beta parameter is provided, the
448  /// initial values of the destination matrix are scaled by beta, and the
449  /// product of the left and right matrices is added to the destination
450  /// matrix. If the alpha parameter is provided, the product of the left
451  /// and right matrices is scaled by alpha.
452  ///
453  /// \param lhs The matrix left of the operator.
454  /// \param rhs The matrix right of the operator.
455  /// \param dst The matrix receiving the product.
456  /// \param alpha The scalar applied to the left matrix.
457  /// \param beta The scalar applied to the right matrix.
458  ///
459  static void gemm(
460  const basic_matrix & lhs,
461  const basic_matrix & rhs,
462  basic_matrix & dst,
463  const value_type alpha = value_type(1),
464  const value_type beta = value_type(0))
465  {
466  assert(lhs.get_width() == rhs.get_height());
467  assert(dst.is_size(lhs.get_height(), rhs.get_width()));
468 
470  CblasRowMajor, // Layout
471  CblasNoTrans, // transa
472  CblasNoTrans, // transb
473  int(lhs.get_height()), // m
474  int(rhs.get_width()), // n
475  int(lhs.get_width()), // k
476  alpha, // alpha
477  lhs.get_data(), // a
478  int(lhs.get_width()), // lda
479  rhs.get_data(), // b
480  int(rhs.get_width()), // ldb
481  beta, // beta
482  dst.get_data(), // c
483  int(dst.get_width())); // ldc
484  }
485 
486  ///
487  /// Computes the solution to the system of linear equations with a
488  /// square coefficient matrix A and multiple right-hand sides.
489  /// Specifically, this method solves for X the system of linear
490  /// equations A*X = B, where A is an n-by-n matrix, the columns of
491  /// matrix B are individual right-hand sides, and the columns of X are
492  /// the corresponding solutions.
493  ///
494  /// \return True if successful; otherwise, false.
495  ///
496  bool gesv()
497  {
498  const auto cx = get_width();
499  const auto cy = get_height();
500 
501  assert(cx > cy);
502 
503  std::vector<int> ipiv;
504  ipiv.resize(cy);
505 
506  return 0 == lapack_type::gesv(
507  lapack_type::row_major, // layout
508  int(cy), // n
509  int(cx - cy), // nrhs
510  get_data(), // a
511  int(cx), // lda
512  ipiv.data(), // ipiv
513  get_data() + cy, // b
514  int(cx)); // ldb
515  }
516 
517  ///
518  /// \return The sum of a column.
519  ///
521  const size_t column) ///< The column to sum.
522  const
523  {
524  assert(column < get_width());
525 
526  auto sum = value_type(0);
527  auto ptr = get_data(0, column);
528  const auto end = ptr + get_length();
529 
530  while (ptr != end)
531  {
532  sum += *ptr;
533  ptr += get_width();
534  }
535 
536  return sum;
537  }
538 
539  ///
540  /// \return A pointer to the first element of the data.
541  ///
542  inline const value_type * get_data() const
543  {
544  assert(!is_empty());
545  return _m.data();
546  }
547 
548  ///
549  /// \return A pointer to the first element of the data.
550  ///
551  inline value_type * get_data()
552  {
553  assert(!is_empty());
554  return _m.data();
555  }
556 
557  ///
558  /// \return A pointer to the specified element of the vector data.
559  ///
560  inline const value_type * get_data(
561  const size_t index) ///< The index of the element.
562  const
563  {
564  assert(is_vector());
565  assert(index < get_length());
566  return _m.data() + index;
567  }
568 
569  ///
570  /// \return A pointer to the specified element of the vector data.
571  ///
573  const size_t index) ///< The index of the element.
574  {
575  assert(index < get_length());
576  return _m.data() + index;
577  }
578 
579  ///
580  /// \return AReturns a pointer to the specified element of the vector data.
581  ///
582  inline const value_type * get_data(
583  const size_t row, ///< The row.
584  const size_t column) ///< The column.
585  const
586  {
587  return _m.data() + get_index(row, column);
588  }
589 
590  ///
591  /// \return A pointer to the specified element of the vector data.
592  ///
594  const size_t row, ///< The row of the element.
595  const size_t column) ///< The column of the element.
596  {
597  return _m.data() + get_index(row, column);
598  }
599 
600  ///
601  /// \return The height of the matrix.
602  ///
603  inline size_t get_height() const
604  {
605  return _cy;
606  }
607 
608  ///
609  /// \return The element index for the specified row and column.
610  ///
611  inline size_t get_index(
612  const size_t row, ///< The row.
613  const size_t column) ///< The column.
614  const
615  {
616  assert(row < _cy);
617  assert(column < _cx);
618  return row * _cx + column;
619  }
620 
621  ///
622  /// \return The length of the matrix.
623  ///
624  size_t get_length() const
625  {
626  return _m.size();
627  }
628 
629  ///
630  /// Returns the maximum element value. This method should not be called
631  /// if the matrix is empty.
632  ///
633  /// \return The maximum element value.
634  ///
636  {
637  assert(!is_empty());
638  return *std::max_element(_m.begin(), _m.end());
639  }
640 
641  ///
642  /// Returns the minimum element value. This method should not be called
643  /// if the matrix is empty.
644  ///
645  /// \return The minimum element value.
646  ///
648  {
649  assert(!is_empty());
650  return *std::min_element(_m.begin(), _m.end());
651  }
652 
653  ///
654  /// Returns the minimum and maximum elements in the matrix. If the
655  /// matrix is empty, this method returns false. Otherwise, this method
656  /// stores the results in the specified arguments and returns true.
657  ///
658  /// \param min The minimum value (output).
659  /// \param max The maximum value (output).
660  ///
661  /// \return True if successful; otherwise, false.
662  ///
663  bool get_min_max(value_type & min, value_type & max) const
664  {
665  if (_m.empty())
666  return false;
667 
668  const auto p = *std::minmax_element(_m.begin(), _m.end());
669  min = p.first;
670  max = p.second;
671  return true;
672  }
673 
674 
675  ///
676  /// Returns the minimum and maximum elements in a column. If the
677  /// matrix is empty, this method returns false. Otherwise, this method
678  /// stores the results in the specified arguments and returns true.
679  ///
680  /// \param column The column index.
681  /// \param min The minimum value (output).
682  /// \param max The maximum value (output).
683  ///
684  /// \return True if successful; otherwise, false.
685  ///
687  const size_t column, ///< The column index.
688  value_type & min, ///< The minimum column element.
689  value_type & max) ///< The maximum column element.
690  const
691  {
692  assert(column < get_width());
693 
694  if (_m.empty())
695  return false;
696 
697  auto ptr = get_data(0, column);
698  auto const end = ptr + get_length();
699 
700  min = max = *ptr;
701 
702  for (ptr += get_width(); ptr < end; ptr += get_width())
703  {
704  const auto v = *ptr;
705  min = std::min(min, v);
706  max = std::max(max, v);
707  }
708 
709  return true;
710  }
711 
712 
713  ///
714  /// \return The sum of a row.
715  ///
717  const size_t row) ///< The row to sum.
718  const
719  {
720  assert(row < get_height());
721 
722  auto sum = value_type(0);
723  auto ptr = get_data(row, 0);
724  const auto end = ptr + get_width();
725 
726  while (ptr < end)
727  sum += *ptr++;
728 
729  return sum;
730  }
731 
732  ///
733  /// \return A string representation of the size of the matrix.
734  ///
735  std::string get_size_str() const
736  {
737  std::ostringstream out;
738  out << '[' << _cy << 'x' << _cx << ']';
739  return out.str();
740  }
741 
742  ///
743  /// \return The sum of all elements in the matrix.
744  ///
746  {
747  return std::accumulate(
748  _m.begin(),
749  _m.end(),
750  value_type(0),
751  std::plus<value_type>());
752  }
753 
754  ///
755  /// \return The height of the matrix.
756  ///
757  inline size_t get_width() const
758  {
759  return _cx;
760  }
761 
762  ///
763  /// \return Thevalue at the specified index in the vector data.
764  ///
766  const size_t index) ///< The index in the vector data.
767  const
768  {
769  assert(index < get_length());
770  return _m[index];
771  }
772 
773  ///
774  /// \return The value at the specified row and column.
775  ///
777  const size_t row, ///< The row index.
778  const size_t column) ///< The column index.
779  const
780  {
781  const auto index = get_index(row, column);
782  return _m[index];
783  }
784 
785  ///
786  /// Computes and stores the inverse of this matrix using the Cholesky
787  /// square root method. Only the lower-left triangle is used to perform
788  /// the calculation. If this method fails, indicating this matrix is
789  /// not positive semidefinite, then the contents of this instance and
790  /// the parameter are left undefined, and the method returns false.
791  /// Otherwise, the lower triangle is copied to the upper triangle,
792  /// the parameter is assigned the log of the determinant of the
793  /// original matrix values, and the method returns true.
794  /// \return True if successful; otherwise, false.
795  ///
796  bool invert(
797  value_type & log_det) ///< The log of the determinant.
798  {
799  assert(is_square());
800 
801  //
802  // Compute the Cholesky square root. If it fails, return false to
803  // indicate this is an unacceptable set of parameters.
804  //
805  if (!potrf_lower())
806  return false;
807 
808  //
809  // Calculate the log of the determinant by summing twice the log of
810  // the diagonal entries.
811  //
812  log_det = value_type(0.0);
813  {
814  const auto end = get_data() + get_length() + _cx;
815  for (auto ptr = get_data(); ptr != end; ptr += _cx + 1)
816  log_det += value_type(2.0) * std::log(*ptr);
817  }
818 
819  //
820  // Compute the inverse. If it fails, return false to indicate
821  // this is an unacceptable set of parameters.
822  //
823  if (!potri_lower())
824  return false;
825 
826  //
827  // Mirror the values from the lower triangle to the upper triangle.
828  //
830  return true;
831  }
832 
833  ///
834  /// Computes and stores the inverse of this matrix using the Cholesky
835  /// square root method. Only the lower-left triangle is used to perform
836  /// the calculation. If this method fails, indicating this matrix is
837  /// not positive semidefinite, then the contents of this instance is
838  /// left undefined, and the method returns false. Otherwise, the lower
839  /// triangle is copied to the upper triangle, and the method returns
840  /// true.
841  /// \return True if successful; otherwise, false.
842  ///
843  inline bool invert()
844  {
845  value_type log_det;
846  return invert(log_det);
847  }
848 
849  ///
850  /// \return True if this is a column vector.
851  ///
852  inline bool is_column_vector() const
853  {
854  return _cx == 1;
855  }
856 
857  ///
858  /// \return True if this is empty.
859  ///
860  inline bool is_empty() const
861  {
862  return _m.empty();
863  }
864 
865  ///
866  /// \return True if this is a row vector.
867  ///
868  inline bool is_row_vector() const
869  {
870  return _cy == 1;
871  }
872 
873  ///
874  /// \return True if this has the same length as a specified vector.
875  ///
876  inline bool is_length(
877  const basic_matrix & other) ///< The other matrix.
878  const
879  {
880  return get_length() == other.get_length();
881  }
882 
883  ///
884  /// \return True if this has the specified length.
885  ///
886  inline bool is_length(
887  const size_t length) ///< The length to compare.
888  const
889  {
890  return get_length() == length;
891  }
892 
893  ///
894  /// \return True if this has the same size as a specified matrix.
895  ///
896  inline bool is_size(
897  const basic_matrix & other) ///< The other matrix.
898  const
899  {
900  return (_cx == other._cx) && (_cy == other._cy);
901  }
902 
903  ///
904  /// \return True if this has the same size as the specified dimensions.
905  ///
906  inline bool is_size(
907  const size_t height, ///< The height of the matrix.
908  const size_t width) ///< The width of the matrix.
909  const
910  {
911  return (_cy == height) && (_cx == width);
912  }
913 
914  ///
915  /// \return True if this matrix is square.
916  ///
917  inline bool is_square() const
918  {
919  return _cx == _cy;
920  }
921 
922  ///
923  /// \return True if this is a vector.
924  ///
925  inline bool is_vector() const
926  {
927  return is_column_vector() || is_row_vector();
928  }
929 
930  ///
931  /// Multiplies a column by a specified value.
932  ///
934  const size_t column, ///< The column to multiply.
935  const value_type value) ///< The value to multiply.
936  {
937  assert(column < get_width());
938 
939  auto dst_ptr = get_data(0, column);
940  const auto dst_end = dst_ptr + get_length();
941 
942  while (dst_ptr != dst_end)
943  {
944  *dst_ptr *= value;
945  dst_ptr += get_width();
946  }
947  }
948 
949  ///
950  /// \return The product of the column and the vector.
951  ///
953  const size_t column, ///< The column to multiply.
954  const basic_matrix & vector) ///< The vector to multiply.
955  {
956  assert(column < get_height());
957  assert(vector.is_vector());
958 
959  auto a_ptr = get_data(0, column);
960  const auto a_end = a_ptr + get_length();
961  auto b_ptr = vector.get_data();
962  auto sum = value_type(0);
963 
964  while (a_ptr != a_end)
965  {
966  sum += *a_ptr * *b_ptr++;
967  a_ptr += get_width();
968  }
969 
970  return sum;
971  }
972 
973  ///
974  /// Multiplies a row by a specified value.
975  ///
977  const size_t row, ///< The row to multiply.
978  const value_type value) ///< The value to multiply.
979  {
980  assert(row < get_height());
981 
982  auto dst_ptr = get_data(row, 0);
983  const auto dst_end = dst_ptr + get_width();
984 
985  while (dst_ptr != dst_end)
986  *dst_ptr++ *= value;
987  }
988 
989  ///
990  /// \return The product of the row and the vector.
991  ///
993  const size_t row, ///< The row to multiply.
994  const basic_matrix & vector) ///< The vector to multiply.
995  const
996  {
997  assert(row < get_height());
998  assert(vector.is_vector());
999 
1000  auto a_ptr = get_data(row, 0);
1001  const auto a_end = a_ptr + get_width();
1002  auto b_ptr = vector.get_data();
1003  auto sum = value_type(0);
1004 
1005  while (a_ptr != a_end)
1006  sum += *a_ptr++ * *b_ptr++;
1007 
1008  return sum;
1009  }
1010 
1011  ///
1012  /// Forms the Cholesky factorization of a symmetric positive-definite
1013  /// matrix: [A = L * L^T] where L is the lower triangular portion of the
1014  /// matrix.
1015  ///
1016  /// \return True if the method is successful.
1017  ///
1019  {
1020  assert(is_square());
1021  assert(!is_empty());
1022 
1023  const auto n = int(get_width());
1024 
1025  return 0 == lapack_type::potrf(
1026  lapack_type::row_major, // layout
1027  'L', // uplo
1028  n, // n
1029  get_data(), // a
1030  n); // lda
1031  }
1032 
1033  ///
1034  /// Computes the inverse inv(A) of a symmetric positive definite matrix
1035  /// A. Before calling this routine, call dpotrf to factorize A.
1036  ///
1037  /// \return True if the method is successful.
1038  ///
1040  {
1041  assert(is_square());
1042  assert(!is_empty());
1043 
1044  const auto n = int(get_width());
1045 
1046  return 0 == lapack_type::potri(
1047  lapack_type::row_major, // layout
1048  'L', // uplo
1049  n, // n
1050  get_data(), // a
1051  n); // lda
1052  }
1053 
1054  ///
1055  /// Reads the matrix values from the specified file.
1056  ///
1057  /// \throw An exception if there is an error reading the file.
1058  ///
1059  inline void read(
1060  const std::string & path) ///< The path to the file.
1061  {
1062  read(path.c_str());
1063  }
1064 
1065  ///
1066  /// Reads the matrix values from the specified file.
1067  ///
1068  /// \throw An exception if there is an error reading the file.
1069  ///
1070  void read(
1071  const char * const path) ///< The path to the file.
1072  {
1073  assert(path != nullptr);
1074 
1075  std::ifstream in (path);
1076  if (!in.good())
1077  throw error() << "failed to open matrix '" << path << "'";
1078 
1079  try
1080  {
1081  read(in);
1082  }
1083  catch (const std::exception & e)
1084  {
1085  throw error()
1086  << "failed to read matrix '"
1087  << path << "': " << e.what();
1088  }
1089  }
1090 
1091  ///
1092  /// Reads the matrix values from the specified stream.
1093  ///
1094  /// \throw An exception if there is an error reading the stream.
1095  ///
1096  void read(
1097  std::istream & in) ///< The input stream.
1098  {
1099  size_t cx, cy;
1100  if (!(in >> cy >> cx))
1101  throw error()
1102  << "failed to parse matrix size";
1103 
1104  if ((cx == 0 && cy != 0) || (cx != 0 && cy == 0))
1105  throw error()
1106  << "invalid matrix size ["
1107  << cy << "x" << cx << "]";
1108 
1109  vector_type m;
1110  m.reserve(cx * cy);
1111 
1112  for (size_t y = 0; y < cy; y++)
1113  {
1114  for (size_t x = 0; x < cx; x++)
1115  {
1116  value_type value;
1117  if (!(in >> value))
1118  throw error()
1119  << "failed to parse matrix value at cell ["
1120  << y+1 << "," << x+1 << "]";
1121 
1122  m.push_back(value);
1123  }
1124  }
1125 
1126  _cx = cx;
1127  _cy = cy;
1128  _m.swap(m);
1129  }
1130 
1131  ///
1132  /// Resizes the matrix to the specified dimensions. The values of the
1133  /// matrix are not reset.
1134  ///
1135  void resize(
1136  const size_t height, ///< The new height.
1137  const size_t width) ///< The new width.
1138  {
1139  assert((height == 0 && width == 0)
1140  || (height != 0 && width != 0));
1141 
1142  if (height == _cy && width == _cx)
1143  return;
1144 
1145  _cx = width;
1146  _cy = height;
1147  _m.resize(width * height);
1148  }
1149 
1150  ///
1151  /// Sets the specified stream to scientific notation with a high
1152  /// precision.
1153  ///
1154  inline static void set_high_precision(
1155  std::ostream & out) ///< The output stream.
1156  {
1157  static const auto precision =
1158  1 + std::numeric_limits<value_type>::digits10;
1159 
1160  out << std::scientific << std::setprecision(precision);
1161  }
1162 
1163  ///
1164  /// Sets a single value of the vector.
1165  ///
1166  inline void set_value(
1167  const size_t index, ///< The index of the value.
1168  const value_type & value) ///< The new value.
1169  {
1170  assert(index < _m.size());
1171  _m[index] = value;
1172  }
1173 
1174  ///
1175  /// Sets a single value of the matrix.
1176  ///
1177  inline void set_value(
1178  const size_t row, ///< The row to assign.
1179  const size_t column, ///< The colum to assign.
1180  const value_type & value) ///< The value to assign.
1181  {
1182  const auto index = get_index(row, column);
1183  _m[index] = value;
1184  }
1185 
1186  ///
1187  /// Sets all values of the matrix to the specified value.
1188  ///
1190  const value_type value) ///< The value to assign.
1191  {
1192  for (auto & m : _m)
1193  m = value;
1194  }
1195 
1196  ///
1197  /// \return A string representation of this instance.
1198  ///
1199  std::string str() const
1200  {
1201  std::ostringstream out;
1202  write(out);
1203  return out.str();
1204  }
1205 
1206  ///
1207  /// Swaps this matrix and another matrix.
1208  ///
1209  inline void swap(
1210  basic_matrix & other) ///< The matrix to swap.
1211  {
1212  std::swap(_cx, other._cx);
1213  std::swap(_cy, other._cy);
1214  _m.swap(other._m);
1215  }
1216 
1217  ///
1218  /// Transposes this matrix.
1219  ///
1220  inline void transpose()
1221  {
1222  create_transpose(*this);
1223  }
1224 
1225  ///
1226  /// Writes this matrix to the specified file.
1227  ///
1228  /// \throw An exception if there is an error writing to the file.
1229  ///
1230  void write(
1231  const std::string & path) ///< The path to the file.
1232  const
1233  {
1234  write(path.c_str());
1235  }
1236 
1237  ///
1238  /// Writes this matrix to the specified file.
1239  ///
1240  /// \throw An exception if there is an error writing to the file.
1241  ///
1242  void write(
1243  const char * const path) ///< The path to the file.
1244  const
1245  {
1246  assert(path != nullptr);
1247 
1248  std::ofstream out (path);
1249 
1250  if (!out.good())
1251  throw error() << "failed to create matrix '" << path << "'";
1252 
1253  write(out);
1254  }
1255 
1256  ///
1257  /// Writes this matrix to the specified output stream.
1258  ///
1259  void write(
1260  std::ostream & out) ///< The output stream.
1261  const
1262  {
1263  value_type const * ptr = _m.data();
1264 
1265  out << _cy << ' ' << _cx << '\n';
1266 
1267  if (_cx == 0)
1268  return;
1269 
1270  for (size_t y = 0; y < _cy; y++)
1271  {
1272  out << *ptr++;
1273 
1274  for (size_t x = 1; x < _cx; x++)
1275  out << '\t' << *ptr++;
1276 
1277  out << '\n';
1278  }
1279  }
1280 
1281  ///
1282  /// \return The element at the specified index.
1283  ///
1284  inline const value_type & operator [] (
1285  const size_t rhs) ///< The element index.
1286  const
1287  {
1288  assert(rhs < get_length());
1289  return _m[rhs];
1290  }
1291 
1292  ///
1293  /// \return The element at the specified index.
1294  ///
1296  const size_t rhs) ///< The element index.
1297  {
1298  assert(rhs < get_length());
1299  return _m[rhs];
1300  }
1301 
1302  ///
1303  /// \return The element at the specified row and column.
1304  ///
1305  inline const value_type & operator () (
1306  const size_t rhs_1, ///< The row.
1307  const size_t rhs_2) ///< The column.
1308  const
1309  {
1310  const auto index = get_index(rhs_1, rhs_2);
1311  return _m[index];
1312  }
1313 
1314  ///
1315  /// \return The element at the specified row and column.
1316  ///
1318  const size_t rhs_1, ///< The row.
1319  const size_t rhs_2) ///< The column.
1320  {
1321  const auto index = get_index(rhs_1, rhs_2);
1322  return _m[index];
1323  }
1324 
1325  ///
1326  /// Adds the specified matrix to this instance.
1327  ///
1328  /// \param rhs The matrix to add.
1329  /// \return This instance.
1330  ///
1332  {
1333  assert(get_width() == rhs.get_width());
1334  assert(get_height() == rhs.get_height());
1335 
1336  auto dst_ptr = get_data();
1337  auto src_ptr = rhs.get_data();
1338  const auto dst_end = dst_ptr + get_length();
1339 
1340  while (dst_ptr != dst_end)
1341  *dst_ptr++ += *src_ptr++;
1342 
1343  return *this;
1344  }
1345 
1346  ///
1347  /// Adds the specified scalar to this instance.
1348  ///
1349  /// \param rhs The scalar to add.
1350  /// \return This instance.
1351  ///
1353  {
1354  auto dst_ptr = get_data();
1355  const auto dst_end = dst_ptr + get_length();
1356 
1357  while (dst_ptr != dst_end)
1358  *dst_ptr++ += rhs;
1359 
1360  return *this;
1361  }
1362 
1363  ///
1364  /// Subtracts the specified matrix from this instance.
1365  ///
1366  /// \param rhs The matrix to subtract.
1367  /// \return This instance.
1368  ///
1370  {
1371  assert(get_width() == rhs.get_width());
1372  assert(get_height() == rhs.get_height());
1373 
1374  auto dst_ptr = get_data();
1375  auto src_ptr = rhs.get_data();
1376  const auto dst_end = dst_ptr + get_length();
1377 
1378  while (dst_ptr != dst_end)
1379  *dst_ptr++ -= *src_ptr++;
1380 
1381  return *this;
1382  }
1383 
1384  ///
1385  /// Subtracts the specified scalar from this instance.
1386  ///
1387  /// \param rhs The scalar to subtract.
1388  /// \return This instance.
1389  ///
1391  {
1392  auto dst_ptr = get_data();
1393  const auto dst_end = dst_ptr + get_length();
1394 
1395  while (dst_ptr != dst_end)
1396  *dst_ptr++ -= rhs;
1397 
1398  return *this;
1399  }
1400 
1401  ///
1402  /// Multiplies this matrix and the specified matrix and assigns the
1403  /// product to this instance.
1404  ///
1405  /// \param rhs The matrix to multiply.
1406  /// \return This instance.
1407  ///
1409  {
1410  basic_matrix dst (get_height(), rhs.get_width());
1411  gemm(*this, rhs, dst);
1412  swap(dst);
1413  return *this;
1414  }
1415 
1416  ///
1417  /// Multiplies this matrix and the specified scalar and assigns the
1418  /// product to this instance.
1419  ///
1420  /// \param rhs The scalar to multiply.
1421  /// \return This instance.
1422  ///
1424  {
1425  auto dst_ptr = get_data();
1426  const auto dst_end = dst_ptr + get_length();
1427 
1428  while (dst_ptr != dst_end)
1429  *dst_ptr++ *= rhs;
1430 
1431  return *this;
1432  }
1433 
1434  ///
1435  /// Divides this matrix by the specified scalar and assigns the
1436  /// quotient to this instance.
1437  ///
1438  /// \param rhs The scalar to multiply.
1439  /// \return This instance.
1440  ///
1442  {
1443  auto dst_ptr = get_data();
1444  const auto dst_end = dst_ptr + get_length();
1445 
1446  while (dst_ptr != dst_end)
1447  *dst_ptr++ /= rhs;
1448 
1449  return *this;
1450  }
1451 
1452  ///
1453  /// Adds a matrix.
1454  ///
1455  /// \param rhs The matrix to add.
1456  /// \return A new matrix.
1457  ///
1458  inline basic_matrix operator + (const basic_matrix & rhs) const
1459  {
1460  return basic_matrix(*this) += rhs;
1461  }
1462 
1463  ///
1464  /// Adds a value.
1465  ///
1466  /// \param rhs The value to add.
1467  /// \return A new matrix.
1468  ///
1469  inline basic_matrix operator + (const value_type rhs) const
1470  {
1471  return basic_matrix(*this) += rhs;
1472  }
1473 
1474  ///
1475  /// Subtracts a matrix.
1476  ///
1477  /// \param rhs The matrix to subtract.
1478  /// \return A new matrix.
1479  ///
1480  inline basic_matrix operator - (const basic_matrix & rhs) const
1481  {
1482  return basic_matrix(*this) -= rhs;
1483  }
1484 
1485  ///
1486  /// Subtracts a value.
1487  ///
1488  /// \param rhs The value to subtract.
1489  /// \return A new matrix.
1490  ///
1491  inline basic_matrix operator - (const value_type rhs) const
1492  {
1493  return basic_matrix(*this) -= rhs;
1494  }
1495 
1496  ///
1497  /// \return The negation of this instance.
1498  ///
1500  {
1501  basic_matrix out (get_height(), get_width());
1502 
1503  auto dst = out.get_data();
1504  auto src = get_data();
1505  const auto end = get_data() + get_length();
1506 
1507  while (src != end)
1508  *dst++ = -(*src++);
1509 
1510  return out;
1511  }
1512 
1513  ///
1514  /// Multiplies a matrix.
1515  ///
1516  /// \param rhs The matrix to multiply.
1517  /// \return A new matrix.
1518  ///
1519  inline basic_matrix operator * (const basic_matrix & rhs) const
1520  {
1521  basic_matrix dst (get_height(), rhs.get_width());
1522  gemm(*this, rhs, dst);
1523  return dst;
1524  }
1525 
1526  ///
1527  /// Multiplies a scalar.
1528  ///
1529  /// \param rhs The scalar to multiply.
1530  /// \return A new matrix.
1531  ///
1532  inline basic_matrix operator * (const value_type rhs) const
1533  {
1534  return basic_matrix(*this) *= rhs;
1535  }
1536 
1537  ///
1538  /// Divides a scalar.
1539  ///
1540  /// \param rhs The scalar to divides.
1541  /// \return A new matrix.
1542  ///
1543  inline basic_matrix operator / (const value_type rhs) const
1544  {
1545  return basic_matrix(*this) /= rhs;
1546  }
1547 
1548  private:
1549  typedef std::vector<value_type> vector_type;
1550 
1551  size_t _cy;
1552  size_t _cx;
1553  vector_type _m;
1554  };
1555 }
1556 
1557 ///
1558 /// Multiplies the specified scalar with the specified matrix.
1559 ///
1560 /// \return A new matrix.
1561 ///
1562 template <typename TValue>
1563 jade::basic_matrix<TValue> operator * (
1564  const TValue lhs, ///< The scalar.
1565  const jade::basic_matrix<TValue> & rhs) ///< The matrix.
1566 {
1567  return rhs * lhs;
1568 }
1569 
1570 ///
1571 /// Writes the specified matrix to the specified output stream.
1572 ///
1573 /// \return The output stream.
1574 ///
1575 template <typename TValue>
1576 std::ostream & operator << (
1577  std::ostream & lhs, ///< The output stream.
1578  const jade::basic_matrix<TValue> & rhs) ///< The matrix.
1579 {
1580  rhs.write(lhs);
1581  return lhs;
1582 }
1583 
1584 ///
1585 /// Reads the specified matrix from the specified output stream.
1586 ///
1587 /// \return The input stream.
1588 ///
1589 template <typename TValue>
1590 std::istream & operator >> (
1591  std::istream & lhs, ///< The input stream.
1592  jade::basic_matrix<TValue> & rhs) ///< The matrix.
1593 {
1594  rhs.read(lhs);
1595  return lhs;
1596 }
1597 
1598 #endif // JADE_MATRIX_HPP__
jade::basic_matrix::operator/
basic_matrix operator/(const value_type rhs) const
Divides a scalar.
Definition: jade.matrix.hpp:1543
jade::basic_matrix::read
void read(const char *const path)
Reads the matrix values from the specified file.
Definition: jade.matrix.hpp:1070
jade::basic_matrix::is_square
bool is_square() const
Definition: jade.matrix.hpp:917
jade::basic_matrix::invert
bool invert()
Computes and stores the inverse of this matrix using the Cholesky square root method....
Definition: jade.matrix.hpp:843
jade::basic_matrix::get_data
const value_type * get_data() const
Definition: jade.matrix.hpp:542
jade::basic_matrix::get_min_max_column
bool get_min_max_column(const size_t column, value_type &min, value_type &max) const
Returns the minimum and maximum elements in a column. If the matrix is empty, this method returns fal...
Definition: jade.matrix.hpp:686
jade::basic_matrix::get_length
size_t get_length() const
Definition: jade.matrix.hpp:624
jade::basic_matrix::all_of
bool all_of(const TPredicate predicate) const
Definition: jade.matrix.hpp:117
jade::basic_matrix::set_value
void set_value(const size_t row, const size_t column, const value_type &value)
Sets a single value of the matrix.
Definition: jade.matrix.hpp:1177
jade::basic_matrix::get_value
value_type get_value(const size_t row, const size_t column) const
Definition: jade.matrix.hpp:776
jade::basic_matrix::create_transpose
void create_transpose(basic_matrix &out) const
Stores the transpose of the matrix.
Definition: jade.matrix.hpp:401
jade::basic_matrix::gesv
bool gesv()
Computes the solution to the system of linear equations with a square coefficient matrix A and multip...
Definition: jade.matrix.hpp:496
jade::basic_matrix::get_index
size_t get_index(const size_t row, const size_t column) const
Definition: jade.matrix.hpp:611
jade::basic_matrix::copy_column
void copy_column(const size_t column, basic_matrix &out) const
Copies a column into the specified vector.
Definition: jade.matrix.hpp:242
jade::basic_matrix::resize
void resize(const size_t height, const size_t width)
Resizes the matrix to the specified dimensions. The values of the matrix are not reset.
Definition: jade.matrix.hpp:1135
jade::basic_matrix::get_width
size_t get_width() const
Definition: jade.matrix.hpp:757
jade::basic_matrix::write
void write(std::ostream &out) const
Writes this matrix to the specified output stream.
Definition: jade.matrix.hpp:1259
jade::basic_matrix::get_data
value_type * get_data()
Definition: jade.matrix.hpp:551
jade::basic_lapack::potrf
static int potrf(layout_type layout, char uplo, int n, value_type *a, int lda)
Computes the Cholesky factorization of a symmetric (Hermitian) positive-definite matrix.
jade::basic_matrix::transpose
void transpose()
Transposes this matrix.
Definition: jade.matrix.hpp:1220
jade::basic_matrix::operator()
const value_type & operator()(const size_t rhs_1, const size_t rhs_2) const
Definition: jade.matrix.hpp:1305
jade::basic_matrix::clamp
void clamp(const value_type min, const value_type max)
Clamps all values in the matrix to the specified range.
Definition: jade.matrix.hpp:138
jade::basic_error::what
virtual const char * what() const
jade::basic_matrix::is_vector
bool is_vector() const
Definition: jade.matrix.hpp:925
jade::basic_matrix::set_high_precision
static void set_high_precision(std::ostream &out)
Sets the specified stream to scientific notation with a high precision.
Definition: jade.matrix.hpp:1154
jade::basic_matrix::swap
void swap(basic_matrix &other)
Swaps this matrix and another matrix.
Definition: jade.matrix.hpp:1209
jade::basic_matrix::basic_matrix
basic_matrix(const initializer_list_type &values)
Initializes a new instance of the class based on the specified values.
Definition: jade.matrix.hpp:98
jade::basic_matrix::multiply_row
value_type multiply_row(const size_t row, const basic_matrix &vector) const
Definition: jade.matrix.hpp:992
jade::basic_matrix::multiply_column
value_type multiply_column(const size_t column, const basic_matrix &vector)
Definition: jade.matrix.hpp:952
jade::basic_matrix::get_data
const value_type * get_data(const size_t row, const size_t column) const
Definition: jade.matrix.hpp:582
jade::basic_matrix::str
std::string str() const
Definition: jade.matrix.hpp:1199
jade::basic_matrix::write
void write(const std::string &path) const
Writes this matrix to the specified file.
Definition: jade.matrix.hpp:1230
jade::basic_matrix::operator-=
basic_matrix & operator-=(const basic_matrix &rhs)
Subtracts the specified matrix from this instance.
Definition: jade.matrix.hpp:1369
jade::basic_matrix::get_size_str
std::string get_size_str() const
Definition: jade.matrix.hpp:735
jade::basic_matrix::is_size
bool is_size(const size_t height, const size_t width) const
Definition: jade.matrix.hpp:906
jade::basic_matrix::contains_nan
bool contains_nan() const
Definition: jade.matrix.hpp:214
jade::basic_matrix::value_type
TValue value_type
The value type.
Definition: jade.matrix.hpp:24
jade::basic_matrix::read
void read(std::istream &in)
Reads the matrix values from the specified stream.
Definition: jade.matrix.hpp:1096
jade::basic_matrix::contains_inf
bool contains_inf() const
Definition: jade.matrix.hpp:199
jade::basic_matrix::get_sum
value_type get_sum() const
Definition: jade.matrix.hpp:745
jade::basic_matrix::multiply_column
void multiply_column(const size_t column, const value_type value)
Multiplies a column by a specified value.
Definition: jade.matrix.hpp:933
jade::basic_matrix::read
void read(const std::string &path)
Reads the matrix values from the specified file.
Definition: jade.matrix.hpp:1059
jade::basic_matrix::is_size
bool is_size(const basic_matrix &other) const
Definition: jade.matrix.hpp:896
jade::basic_matrix::operator[]
const value_type & operator[](const size_t rhs) const
Definition: jade.matrix.hpp:1284
jade::basic_blas
A template for a class providing access to BLAS.
Definition: jade.blas.hpp:19
jade::basic_matrix::blas_type
basic_blas< value_type > blas_type
The BLAS type.
Definition: jade.matrix.hpp:27
jade::basic_matrix::potrf_lower
bool potrf_lower()
Forms the Cholesky factorization of a symmetric positive-definite matrix: [A = L * L^T] where L is th...
Definition: jade.matrix.hpp:1018
jade::basic_matrix::get_data
value_type * get_data(const size_t row, const size_t column)
Definition: jade.matrix.hpp:593
jade::basic_lapack::potri
static int potri(layout_type layout, char uplo, int n, value_type *a, int lda)
Computes the inverse of a symmetric (Hermitian) positive-definite matrix using the Cholesky factoriza...
jade::basic_matrix::lapack_type
basic_lapack< value_type > lapack_type
The LAPACK type.
Definition: jade.matrix.hpp:30
jade::basic_matrix::get_value
value_type get_value(const size_t index) const
Definition: jade.matrix.hpp:765
jade::basic_matrix::is_empty
bool is_empty() const
Definition: jade.matrix.hpp:860
jade::basic_matrix::operator*
basic_matrix operator*(const basic_matrix &rhs) const
Multiplies a matrix.
Definition: jade.matrix.hpp:1519
jade::basic_matrix::basic_matrix
basic_matrix(const std::string &path)
Initializes a new instance of the class based on values from the specified file.
Definition: jade.matrix.hpp:75
jade::basic_matrix::write
void write(const char *const path) const
Writes this matrix to the specified file.
Definition: jade.matrix.hpp:1242
jade::basic_matrix::basic_matrix
basic_matrix(const size_t cy, const size_t cx)
Initializes a new instance of the class based on the specified width and height.
Definition: jade.matrix.hpp:49
jade::basic_matrix::basic_matrix
basic_matrix()
Initializes a new instance of the class with no size.
Definition: jade.matrix.hpp:40
jade::basic_matrix::get_height
size_t get_height() const
Definition: jade.matrix.hpp:603
jade::basic_lapack::gesv
static int gesv(layout_type layout, int n, int nrhs, value_type *a, int lda, int *ipiv, value_type *b, int ldb)
Computes the solution to the system of linear equations with a square coefficient matrix A and multip...
jade::basic_matrix::basic_matrix
basic_matrix(std::istream &in)
Initializes a new instance of the class based on values from the specified input stream.
Definition: jade.matrix.hpp:85
jade::basic_matrix::get_min_max
bool get_min_max(value_type &min, value_type &max) const
Returns the minimum and maximum elements in the matrix. If the matrix is empty, this method returns f...
Definition: jade.matrix.hpp:663
jade::basic_matrix::invert
bool invert(value_type &log_det)
Computes and stores the inverse of this matrix using the Cholesky square root method....
Definition: jade.matrix.hpp:796
jade::basic_matrix::operator+=
basic_matrix & operator+=(const basic_matrix &rhs)
Adds the specified matrix to this instance.
Definition: jade.matrix.hpp:1331
jade::basic_matrix::get_max_value
value_type get_max_value() const
Returns the maximum element value. This method should not be called if the matrix is empty.
Definition: jade.matrix.hpp:635
jade::basic_matrix::clamp_column
void clamp_column(const size_t column, const value_type min, const value_type max)
Clamps all values in a column to the specified range.
Definition: jade.matrix.hpp:157
jade::basic_matrix::set_value
void set_value(const size_t index, const value_type &value)
Sets a single value of the vector.
Definition: jade.matrix.hpp:1166
jade::basic_matrix::is_length
bool is_length(const size_t length) const
Definition: jade.matrix.hpp:886
jade::basic_matrix::copy_upper_to_lower
void copy_upper_to_lower()
Copies the elements in the upper triangle of the square matrix to the corresponding lower triangle el...
Definition: jade.matrix.hpp:345
jade::basic_matrix::get_data
value_type * get_data(const size_t index)
Definition: jade.matrix.hpp:572
jade::basic_matrix::multiply_row
void multiply_row(const size_t row, const value_type value)
Multiplies a row by a specified value.
Definition: jade.matrix.hpp:976
jade::basic_matrix::get_column_sum
value_type get_column_sum(const size_t column) const
Definition: jade.matrix.hpp:520
jade::basic_matrix::initializer_list_type
std::initializer_list< std::initializer_list< value_type > > initializer_list_type
The initializer list type.
Definition: jade.matrix.hpp:35
jade::basic_matrix::is_length
bool is_length(const basic_matrix &other) const
Definition: jade.matrix.hpp:876
jade::basic_matrix::operator/=
basic_matrix & operator/=(const value_type rhs)
Divides this matrix by the specified scalar and assigns the quotient to this instance.
Definition: jade.matrix.hpp:1441
jade::basic_matrix::is_row_vector
bool is_row_vector() const
Definition: jade.matrix.hpp:868
jade::basic_matrix::operator-
basic_matrix operator-() const
Definition: jade.matrix.hpp:1499
jade::basic_matrix::operator*=
basic_matrix & operator*=(const basic_matrix &rhs)
Multiplies this matrix and the specified matrix and assigns the product to this instance.
Definition: jade.matrix.hpp:1408
jade::basic_matrix::basic_matrix
basic_matrix(const char *const path)
Initializes a new instance of the class based on values from the specified file.
Definition: jade.matrix.hpp:64
jade::basic_matrix::get_row_sum
value_type get_row_sum(const size_t row) const
Definition: jade.matrix.hpp:716
jade::basic_matrix::is_column_vector
bool is_column_vector() const
Definition: jade.matrix.hpp:852
jade::basic_matrix::copy_lower_to_upper
void copy_lower_to_upper()
Copies the elements in the lower triangle of the square matrix to the corresponding upper triangle el...
Definition: jade.matrix.hpp:266
jade::basic_error
A template for a class representing an exception thrown from this namespace.
Definition: jade.error.hpp:20
jade::basic_matrix::set_values
void set_values(const value_type value)
Sets all values of the matrix to the specified value.
Definition: jade.matrix.hpp:1189
jade::basic_blas::gemm
static void gemm(const layout_type Order, const transpose_type TransA, const transpose_type TransB, const int M, const int N, const int K, const value_type alpha, const value_type *A, const int lda, const value_type *B, const int ldb, const value_type beta, value_type *C, const int ldc)
Computes a matrix-matrix product with general matrices.
jade::basic_matrix::any_of
bool any_of(const TPredicate predicate) const
Definition: jade.matrix.hpp:128
jade::basic_matrix::operator+
basic_matrix operator+(const basic_matrix &rhs) const
Adds a matrix.
Definition: jade.matrix.hpp:1458
jade::basic_matrix::copy_column
basic_matrix copy_column(const size_t column) const
Definition: jade.matrix.hpp:230
jade::basic_matrix::get_data
const value_type * get_data(const size_t index) const
Definition: jade.matrix.hpp:560
jade::basic_matrix::get_min_value
value_type get_min_value() const
Returns the minimum element value. This method should not be called if the matrix is empty.
Definition: jade.matrix.hpp:647
jade::basic_matrix
A template class for a matrix.
Definition: jade.matrix.hpp:21
jade::basic_matrix::copy_row
basic_matrix copy_row(const size_t row) const
Definition: jade.matrix.hpp:312
jade::basic_matrix::potri_lower
bool potri_lower()
Computes the inverse inv(A) of a symmetric positive definite matrix A. Before calling this routine,...
Definition: jade.matrix.hpp:1039
jade::basic_lapack
A template for a class providing access to LAPACK.
Definition: jade.lapack.hpp:19
jade::basic_matrix::clamp_row
void clamp_row(const size_t row, const value_type min, const value_type max)
Clamps all values in a row to the specified range.
Definition: jade.matrix.hpp:178
jade::basic_matrix::gemm
static void gemm(const basic_matrix &lhs, const basic_matrix &rhs, basic_matrix &dst, const value_type alpha=value_type(1), const value_type beta=value_type(0))
Multiplies a left matrix by a right matrix and stores the result into a destination matrix....
Definition: jade.matrix.hpp:459
jade::basic_matrix::create_transpose
basic_matrix create_transpose() const
Definition: jade.matrix.hpp:391
jade::basic_lapack::row_major
@ row_major
A row-major order.
Definition: jade.lapack.hpp:37
jade::basic_matrix::copy_row
void copy_row(const size_t row, basic_matrix &out) const
Copies a row into the specified matrix.
Definition: jade.matrix.hpp:324