GteMatrix.h
Go to the documentation of this file.
1 // David Eberly, Geometric Tools, Redmond WA 98052
2 // Copyright (c) 1998-2017
3 // Distributed under the Boost Software License, Version 1.0.
4 // http://www.boost.org/LICENSE_1_0.txt
5 // http://www.geometrictools.com/License/Boost/LICENSE_1_0.txt
6 // File Version: 3.0.1 (2016/08/07)
7 
8 #pragma once
9 
10 #include <Mathematics/GteVector.h>
11 #include <Mathematics/GteGaussianElimination.h>
12 
13 namespace gte
14 {
15 
16 template <int NumRows, int NumCols, typename Real>
17 class Matrix
18 {
19 public:
20  // The table is initialized to zero.
21  Matrix();
22 
23  // The table is fully initialized by the inputs. The 'values' must be
24  // specified in row-major order, regardless of the active storage scheme
25  // (GTE_USE_ROW_MAJOR or GTE_USE_COL_MAJOR).
26  Matrix(std::array<Real, NumRows*NumCols> const& values);
27 
28  // At most NumRows*NumCols are copied from the initializer list, setting
29  // any remaining elements to zero. The 'values' must be specified in
30  // row-major order, regardless of the active storage scheme
31  // (GTE_USE_ROW_MAJOR or GTE_USE_COL_MAJOR). Create the zero matrix using
32  // the syntax
33  // Matrix<NumRows,NumCols,Real> zero{(Real)0};
34  // WARNING: The C++ 11 specification states that
35  // Matrix<NumRows,NumCols,Real> zero{};
36  // will lead to a call of the default constructor, not the initializer
37  // constructor!
38  Matrix(std::initializer_list<Real> values);
39 
40  // For 0 <= r < NumRows and 0 <= c < NumCols, element (r,c) is 1 and all
41  // others are 0. If either of r or c is invalid, the zero matrix is
42  // created. This is a convenience for creating the standard Euclidean
43  // basis matrices; see also MakeUnit(int,int) and Unit(int,int).
44  Matrix(int r, int c);
45 
46  // The copy constructor, destructor, and assignment operator are generated
47  // by the compiler.
48 
49  // Member access for which the storage representation is transparent. The
50  // matrix entry in row r and column c is A(r,c). The first operator()
51  // returns a const reference rather than a Real value. This supports
52  // writing via standard file operations that require a const pointer to
53  // data.
54  inline Real const& operator()(int r, int c) const;
55  inline Real& operator()(int r, int c);
56 
57  // Member access by rows or by columns.
58  void SetRow(int r, Vector<NumCols,Real> const& vec);
59  void SetCol(int c, Vector<NumRows,Real> const& vec);
60  Vector<NumCols,Real> GetRow(int r) const;
61  Vector<NumRows,Real> GetCol(int c) const;
62 
63  // Member access by 1-dimensional index. NOTE: These accessors are
64  // useful for the manipulation of matrix entries when it does not
65  // matter whether storage is row-major or column-major. Do not use
66  // constructs such as M[c+NumCols*r] or M[r+NumRows*c] that expose the
67  // storage convention.
68  inline Real const& operator[](int i) const;
69  inline Real& operator[](int i);
70 
71  // Comparisons for sorted containers and geometric ordering.
72  inline bool operator==(Matrix const& mat) const;
73  inline bool operator!=(Matrix const& mat) const;
74  inline bool operator< (Matrix const& mat) const;
75  inline bool operator<=(Matrix const& mat) const;
76  inline bool operator> (Matrix const& mat) const;
77  inline bool operator>=(Matrix const& mat) const;
78 
79  // Special matrices.
80  void MakeZero(); // All components are 0.
81  void MakeUnit(int r, int c); // Component (r,c) is 1, all others zero.
82  void MakeIdentity(); // Diagonal entries 1, others 0, even when nonsquare
83  static Matrix Zero();
84  static Matrix Unit(int r, int c);
85  static Matrix Identity();
86 
87 protected:
88  class Table
89  {
90  public:
91  // Storage-order-independent element access as 2D array.
92  inline Real const& operator()(int r, int c) const;
93  inline Real& operator()(int r, int c);
94 
95  // Element access as 1D array. Use this internally only when
96  // the 2D storage order is not relevant.
97  inline Real const& operator[](int i) const;
98  inline Real& operator[](int i);
99 
100 #if defined(GTE_USE_ROW_MAJOR)
101  std::array<std::array<Real,NumCols>,NumRows> mStorage;
102 #else
103  std::array<std::array<Real,NumRows>,NumCols> mStorage;
104 #endif
105  };
106 
107  Table mTable;
108 };
109 
110 // Unary operations.
111 template <int NumRows, int NumCols, typename Real>
112 Matrix<NumRows,NumCols,Real>
113 operator+(Matrix<NumRows,NumCols,Real> const& M);
114 
115 template <int NumRows, int NumCols, typename Real>
116 Matrix<NumRows,NumCols,Real>
117 operator-(Matrix<NumRows,NumCols,Real> const& M);
118 
119 // Linear-algebraic operations.
120 template <int NumRows, int NumCols, typename Real>
121 Matrix<NumRows,NumCols,Real>
122 operator+(
123  Matrix<NumRows,NumCols,Real> const& M0,
124  Matrix<NumRows,NumCols,Real> const& M1);
125 
126 template <int NumRows, int NumCols, typename Real>
127 Matrix<NumRows,NumCols,Real>
128 operator-(
129  Matrix<NumRows,NumCols,Real> const& M0,
130  Matrix<NumRows,NumCols,Real> const& M1);
131 
132 template <int NumRows, int NumCols, typename Real>
133 Matrix<NumRows,NumCols,Real>
134 operator*(Matrix<NumRows,NumCols,Real> const& M, Real scalar);
135 
136 template <int NumRows, int NumCols, typename Real>
137 Matrix<NumRows,NumCols,Real>
138 operator*(Real scalar, Matrix<NumRows,NumCols,Real> const& M);
139 
140 template <int NumRows, int NumCols, typename Real>
141 Matrix<NumRows,NumCols,Real>
142 operator/(Matrix<NumRows,NumCols,Real> const& M, Real scalar);
143 
144 template <int NumRows, int NumCols, typename Real>
145 Matrix<NumRows,NumCols,Real>&
146 operator+=(
147  Matrix<NumRows,NumCols,Real>& M0,
148  Matrix<NumRows,NumCols,Real> const& M1);
149 
150 template <int NumRows, int NumCols, typename Real>
151 Matrix<NumRows,NumCols,Real>&
152 operator-=(
153  Matrix<NumRows,NumCols,Real>& M0,
154  Matrix<NumRows,NumCols,Real> const& M1);
155 
156 template <int NumRows, int NumCols, typename Real>
157 Matrix<NumRows,NumCols,Real>&
158 operator*=(Matrix<NumRows,NumCols,Real>& M, Real scalar);
159 
160 template <int NumRows, int NumCols, typename Real>
161 Matrix<NumRows,NumCols,Real>&
162 operator/=(Matrix<NumRows,NumCols,Real>& M, Real scalar);
163 
164 // Geometric operations.
165 template <int NumRows, int NumCols, typename Real>
166 Real L1Norm(Matrix<NumRows,NumCols,Real> const& M);
167 
168 template <int NumRows, int NumCols, typename Real>
169 Real L2Norm(Matrix<NumRows,NumCols,Real> const& M);
170 
171 template <int NumRows, int NumCols, typename Real>
172 Real LInfinityNorm(Matrix<NumRows,NumCols,Real> const& M);
173 
174 template <int N, typename Real>
175 Matrix<N,N,Real> Inverse(Matrix<N,N,Real> const& M,
176  bool* reportInvertibility = nullptr);
177 
178 template <int N, typename Real>
179 Real Determinant(Matrix<N, N, Real> const& M);
180 
181 // M^T
182 template <int NumRows, int NumCols, typename Real>
183 Matrix<NumCols,NumRows,Real>
184 Transpose(Matrix<NumRows,NumCols,Real> const& M);
185 
186 // M*V
187 template <int NumRows, int NumCols, typename Real>
188 Vector<NumRows,Real>
189 operator*(
190  Matrix<NumRows,NumCols,Real> const& M,
191  Vector<NumCols,Real> const& V);
192 
193 // V^T*M
194 template <int NumRows, int NumCols, typename Real>
195 Vector<NumCols,Real>
196 operator*(
197  Vector<NumRows,Real> const& V,
198  Matrix<NumRows,NumCols,Real> const& M);
199 
200 // A*B
201 template <int NumRows, int NumCols, int NumCommon, typename Real>
202 Matrix<NumRows,NumCols,Real>
203 operator*(
204  Matrix<NumRows,NumCommon,Real> const& A,
205  Matrix<NumCommon,NumCols,Real> const& B);
206 
207 template <int NumRows, int NumCols, int NumCommon, typename Real>
208 Matrix<NumRows,NumCols,Real>
209 MultiplyAB(
210  Matrix<NumRows,NumCommon,Real> const& A,
211  Matrix<NumCommon,NumCols,Real> const& B);
212 
213 // A*B^T
214 template <int NumRows, int NumCols, int NumCommon, typename Real>
215 Matrix<NumRows,NumCols,Real>
216 MultiplyABT(
217  Matrix<NumRows,NumCommon,Real> const& A,
218  Matrix<NumCols,NumCommon,Real> const& B);
219 
220 // A^T*B
221 template <int NumRows, int NumCols, int NumCommon, typename Real>
222 Matrix<NumRows,NumCols,Real>
223 MultiplyATB(
224  Matrix<NumCommon,NumRows,Real> const& A,
225  Matrix<NumCommon,NumCols,Real> const& B);
226 
227 // A^T*B^T
228 template <int NumRows, int NumCols, int NumCommon, typename Real>
229 Matrix<NumRows,NumCols,Real>
230 MultiplyATBT(
231  Matrix<NumCommon,NumRows,Real> const& A,
232  Matrix<NumCols,NumCommon,Real> const& B);
233 
234 // M*D, D is diagonal NumCols-by-NumCols
235 template <int NumRows, int NumCols, typename Real>
236 Matrix<NumRows,NumCols,Real>
237 MultiplyMD(
238  Matrix<NumRows,NumCols,Real> const& M,
239  Vector<NumCols,Real> const& D);
240 
241 // D*M, D is diagonal NumRows-by-NumRows
242 template <int NumRows, int NumCols, typename Real>
243 Matrix<NumRows,NumCols,Real>
244 MultiplyDM(
245  Vector<NumRows,Real> const& D,
246  Matrix<NumRows,NumCols,Real> const& M);
247 
248 // U*V^T, U is NumRows-by-1, V is Num-Cols-by-1, result is NumRows-by-NumCols.
249 template <int NumRows, int NumCols, typename Real>
250 Matrix<NumRows,NumCols,Real>
251 OuterProduct(Vector<NumRows, Real> const& U, Vector<NumCols, Real> const& V);
252 
253 // Initialization to a diagonal matrix whose diagonal entries are the
254 // components of D.
255 template <int N, typename Real>
256 void MakeDiagonal(Vector<N, Real> const& D, Matrix<N, N, Real>& M);
257 
258 // Create an (N+1)-by-(N+1) matrix H by setting the upper N-by-N block to the
259 // input N-by-N matrix and all other entries to 0 except for the last row
260 // and last column entry which is set to 1.
261 template <int N, typename Real>
262 Matrix<N + 1, N+1, Real> HLift(Matrix<N, N, Real> const& M);
263 
264 // Extract the upper (N-1)-by-(N-1) block of the input N-by-N matrix.
265 template <int N, typename Real>
266 Matrix<N - 1, N - 1, Real> HProject(Matrix<N, N, Real> const& M);
267 
268 
269 template <int NumRows, int NumCols, typename Real>
270 Matrix<NumRows, NumCols, Real>::Matrix()
271 {
272  MakeZero();
273 }
274 
275 template <int NumRows, int NumCols, typename Real>
276 Matrix<NumRows, NumCols, Real>::Matrix(
277  std::array<Real, NumRows*NumCols> const& values)
278 {
279  for (int r = 0, i = 0; r < NumRows; ++r)
280  {
281  for (int c = 0; c < NumCols; ++c, ++i)
282  {
283  mTable(r, c) = values[i];
284  }
285  }
286 }
287 
288 template <int NumRows, int NumCols, typename Real>
289 Matrix<NumRows, NumCols, Real>::Matrix(std::initializer_list<Real> values)
290 {
291  int const numValues = static_cast<int>(values.size());
292  auto iter = values.begin();
293  int r, c, i;
294  for (r = 0, i = 0; r < NumRows; ++r)
295  {
296  for (c = 0; c < NumCols; ++c, ++i)
297  {
298  if (i < numValues)
299  {
300  mTable(r, c) = *iter++;
301  }
302  else
303  {
304  break;
305  }
306  }
307 
308  if (c < NumCols)
309  {
310  // Fill in the remaining columns of the current row with zeros.
311  for (; c < NumCols; ++c)
312  {
313  mTable(r, c) = (Real)0;
314  }
315  ++r;
316  break;
317  }
318  }
319 
320  if (r < NumRows)
321  {
322  // Fill in the remain rows with zeros.
323  for (; r < NumRows; ++r)
324  {
325  for (c = 0; c < NumCols; ++c)
326  {
327  mTable(r, c) = (Real)0;
328  }
329  }
330  }
331 }
332 
333 template <int NumRows, int NumCols, typename Real>
334 Matrix<NumRows, NumCols, Real>::Matrix(int r, int c)
335 {
336  MakeUnit(r, c);
337 }
338 
339 template <int NumRows, int NumCols, typename Real> inline
340 Real const& Matrix<NumRows, NumCols, Real>::operator()(int r, int c) const
341 {
342  return mTable(r, c);
343 }
344 
345 template <int NumRows, int NumCols, typename Real> inline
346 Real& Matrix<NumRows, NumCols, Real>::operator()(int r, int c)
347 {
348  return mTable(r, c);
349 }
350 
351 template <int NumRows, int NumCols, typename Real>
352 void Matrix<NumRows, NumCols, Real>::SetRow(int r,
353  Vector<NumCols, Real> const& vec)
354 {
355  for (int c = 0; c < NumCols; ++c)
356  {
357  mTable(r, c) = vec[c];
358  }
359 }
360 
361 template <int NumRows, int NumCols, typename Real>
362 void Matrix<NumRows, NumCols, Real>::SetCol(int c,
363  Vector<NumRows, Real> const& vec)
364 {
365  for (int r = 0; r < NumRows; ++r)
366  {
367  mTable(r, c) = vec[r];
368  }
369 }
370 
371 template <int NumRows, int NumCols, typename Real>
372 Vector<NumCols, Real> Matrix<NumRows, NumCols, Real>::GetRow(int r) const
373 {
374  Vector<NumCols, Real> vec;
375  for (int c = 0; c < NumCols; ++c)
376  {
377  vec[c] = mTable(r, c);
378  }
379  return vec;
380 }
381 
382 template <int NumRows, int NumCols, typename Real>
383 Vector<NumRows, Real> Matrix<NumRows, NumCols, Real>::GetCol(int c) const
384 {
385  Vector<NumRows, Real> vec;
386  for (int r = 0; r < NumRows; ++r)
387  {
388  vec[r] = mTable(r, c);
389  }
390  return vec;
391 }
392 
393 template <int NumRows, int NumCols, typename Real> inline
394 Real const& Matrix<NumRows, NumCols, Real>::operator[](int i) const
395 {
396  return mTable[i];
397 }
398 
399 template <int NumRows, int NumCols, typename Real> inline
400 Real& Matrix<NumRows, NumCols, Real>::operator[](int i)
401 {
402  return mTable[i];
403 }
404 
405 template <int NumRows, int NumCols, typename Real> inline
406 bool Matrix<NumRows, NumCols, Real>::operator==(Matrix const& mat) const
407 {
408  return mTable.mStorage == mat.mTable.mStorage;
409 }
410 
411 template <int NumRows, int NumCols, typename Real> inline
412 bool Matrix<NumRows, NumCols, Real>::operator!=(Matrix const& mat) const
413 {
414  return mTable.mStorage != mat.mTable.mStorage;
415 }
416 
417 template <int NumRows, int NumCols, typename Real> inline
418 bool Matrix<NumRows, NumCols, Real>::operator<(Matrix const& mat) const
419 {
420  return mTable.mStorage < mat.mTable.mStorage;
421 }
422 
423 template <int NumRows, int NumCols, typename Real> inline
424 bool Matrix<NumRows, NumCols, Real>::operator<=(Matrix const& mat) const
425 {
426  return mTable.mStorage <= mat.mTable.mStorage;
427 }
428 
429 template <int NumRows, int NumCols, typename Real> inline
430 bool Matrix<NumRows, NumCols, Real>::operator>(Matrix const& mat) const
431 {
432  return mTable.mStorage > mat.mTable.mStorage;
433 }
434 
435 template <int NumRows, int NumCols, typename Real> inline
436 bool Matrix<NumRows, NumCols, Real>::operator>=(Matrix const& mat) const
437 {
438  return mTable.mStorage >= mat.mTable.mStorage;
439 }
440 
441 template <int NumRows, int NumCols, typename Real>
442 void Matrix<NumRows, NumCols, Real>::MakeZero()
443 {
444  Real const zero = (Real)0;
445  for (int i = 0; i < NumRows * NumCols; ++i)
446  {
447  mTable[i] = zero;
448  }
449 }
450 
451 template <int NumRows, int NumCols, typename Real>
452 void Matrix<NumRows, NumCols, Real>::MakeUnit(int r, int c)
453 {
454  MakeZero();
455  if (0 <= r && r < NumRows && 0 <= c && c < NumCols)
456  {
457  mTable(r, c) = (Real)1;
458  }
459 }
460 
461 template <int NumRows, int NumCols, typename Real>
462 void Matrix<NumRows, NumCols, Real>::MakeIdentity()
463 {
464  MakeZero();
465  int const numDiagonal = (NumRows <= NumCols ? NumRows : NumCols);
466  for (int i = 0; i < numDiagonal; ++i)
467  {
468  mTable(i, i) = (Real)1;
469  }
470 }
471 
472 template <int NumRows, int NumCols, typename Real>
473 Matrix<NumRows, NumCols, Real> Matrix<NumRows, NumCols, Real>::Zero()
474 {
475  Matrix M;
476  M.MakeZero();
477  return M;
478 }
479 
480 template <int NumRows, int NumCols, typename Real>
481 Matrix<NumRows, NumCols, Real> Matrix<NumRows, NumCols, Real>::Unit(int r,
482  int c)
483 {
484  Matrix M;
485  M.MakeUnit(r, c);
486  return M;
487 }
488 
489 template <int NumRows, int NumCols, typename Real>
490 Matrix<NumRows, NumCols, Real> Matrix<NumRows, NumCols, Real>::Identity()
491 {
492  Matrix M;
493  M.MakeIdentity();
494  return M;
495 }
496 
497 
498 
499 template <int NumRows, int NumCols, typename Real>
500 Matrix<NumRows, NumCols, Real>
501 operator+(Matrix<NumRows, NumCols, Real> const& M)
502 {
503  return M;
504 }
505 
506 template <int NumRows, int NumCols, typename Real>
507 Matrix<NumRows, NumCols, Real>
508 operator-(Matrix<NumRows, NumCols, Real> const& M)
509 {
510  Matrix<NumRows, NumCols, Real> result;
511  for (int i = 0; i < NumRows*NumCols; ++i)
512  {
513  result[i] = -M[i];
514  }
515  return result;
516 }
517 
518 template <int NumRows, int NumCols, typename Real>
519 Matrix<NumRows, NumCols, Real>
520 operator+(
521  Matrix<NumRows, NumCols, Real> const& M0,
522  Matrix<NumRows, NumCols, Real> const& M1)
523 {
524  Matrix<NumRows, NumCols, Real> result = M0;
525  return result += M1;
526 }
527 
528 template <int NumRows, int NumCols, typename Real>
529 Matrix<NumRows, NumCols, Real>
530 operator-(
531  Matrix<NumRows, NumCols, Real> const& M0,
532  Matrix<NumRows, NumCols, Real> const& M1)
533 {
534  Matrix<NumRows, NumCols, Real> result = M0;
535  return result -= M1;
536 }
537 
538 template <int NumRows, int NumCols, typename Real>
539 Matrix<NumRows, NumCols, Real>
540 operator*(Matrix<NumRows, NumCols, Real> const& M, Real scalar)
541 {
542  Matrix<NumRows, NumCols, Real> result = M;
543  return result *= scalar;
544 }
545 
546 template <int NumRows, int NumCols, typename Real>
547 Matrix<NumRows, NumCols, Real>
548 operator*(Real scalar, Matrix<NumRows, NumCols, Real> const& M)
549 {
550  Matrix<NumRows, NumCols, Real> result = M;
551  return result *= scalar;
552 }
553 
554 template <int NumRows, int NumCols, typename Real>
555 Matrix<NumRows, NumCols, Real>
556 operator/(Matrix<NumRows, NumCols, Real> const& M, Real scalar)
557 {
558  Matrix<NumRows, NumCols, Real> result = M;
559  return result /= scalar;
560 }
561 
562 template <int NumRows, int NumCols, typename Real>
563 Matrix<NumRows, NumCols, Real>&
564 operator+=(
565  Matrix<NumRows, NumCols, Real>& M0,
566  Matrix<NumRows, NumCols, Real> const& M1)
567 {
568  for (int i = 0; i < NumRows*NumCols; ++i)
569  {
570  M0[i] += M1[i];
571  }
572  return M0;
573 }
574 
575 template <int NumRows, int NumCols, typename Real>
576 Matrix<NumRows, NumCols, Real>&
577 operator-=(
578  Matrix<NumRows, NumCols, Real>& M0,
579  Matrix<NumRows, NumCols, Real> const& M1)
580 {
581  for (int i = 0; i < NumRows*NumCols; ++i)
582  {
583  M0[i] -= M1[i];
584  }
585  return M0;
586 }
587 
588 template <int NumRows, int NumCols, typename Real>
589 Matrix<NumRows, NumCols, Real>&
590 operator*=(Matrix<NumRows, NumCols, Real>& M, Real scalar)
591 {
592  for (int i = 0; i < NumRows*NumCols; ++i)
593  {
594  M[i] *= scalar;
595  }
596  return M;
597 }
598 
599 template <int NumRows, int NumCols, typename Real>
600 Matrix<NumRows, NumCols, Real>&
601 operator/=(Matrix<NumRows, NumCols, Real>& M, Real scalar)
602 {
603  if (scalar != (Real)0)
604  {
605  Real invScalar = ((Real)1) / scalar;
606  for (int i = 0; i < NumRows*NumCols; ++i)
607  {
608  M[i] *= invScalar;
609  }
610  }
611  else
612  {
613  for (int i = 0; i < NumRows*NumCols; ++i)
614  {
615  M[i] = (Real)0;
616  }
617  }
618  return M;
619 }
620 
621 template <int NumRows, int NumCols, typename Real>
622 Real L1Norm(Matrix<NumRows, NumCols, Real> const& M)
623 {
624  Real sum = std::abs(M[0]);
625  for (int i = 1; i < NumRows*NumCols; ++i)
626  {
627  sum += std::abs(M[i]);
628  }
629  return sum;
630 }
631 
632 template <int NumRows, int NumCols, typename Real>
633 Real L2Norm(Matrix<NumRows, NumCols, Real> const& M)
634 {
635  Real sum = M[0] * M[0];
636  for (int i = 1; i < NumRows*NumCols; ++i)
637  {
638  sum += M[i] * M[i];
639  }
640  return sqrt(sum);
641 }
642 
643 template <int NumRows, int NumCols, typename Real>
644 Real LInfinityNorm(Matrix<NumRows, NumCols, Real> const& M)
645 {
646  Real maxAbsElement = M[0];
647  for (int i = 1; i < NumRows*NumCols; ++i)
648  {
649  Real absElement = std::abs(M[i]);
650  if (absElement > maxAbsElement)
651  {
652  maxAbsElement = absElement;
653  }
654  }
655  return maxAbsElement;
656 }
657 
658 template <int N, typename Real>
659 Matrix<N, N, Real> Inverse(Matrix<N, N, Real> const& M,
660  bool* reportInvertibility)
661 {
662  Matrix<N, N, Real> invM;
663  Real determinant;
664  bool invertible = GaussianElimination<Real>()(N, &M[0], &invM[0],
665  determinant, nullptr, nullptr, nullptr, 0, nullptr);
666  if (reportInvertibility)
667  {
668  *reportInvertibility = invertible;
669  }
670  return invM;
671 }
672 
673 template <int N, typename Real>
674 Real Determinant(Matrix<N, N, Real> const& M)
675 {
676  Real determinant;
677  GaussianElimination<Real>()(N, &M[0], nullptr, determinant, nullptr,
678  nullptr, nullptr, 0, nullptr);
679  return determinant;
680 }
681 
682 template <int NumRows, int NumCols, typename Real>
683 Matrix<NumCols, NumRows, Real>
684 Transpose(Matrix<NumRows, NumCols, Real> const& M)
685 {
686  Matrix<NumCols, NumRows, Real> result;
687  for (int r = 0; r < NumRows; ++r)
688  {
689  for (int c = 0; c < NumCols; ++c)
690  {
691  result(c, r) = M(r, c);
692  }
693  }
694  return result;
695 }
696 
697 template <int NumRows, int NumCols, typename Real>
698 Vector<NumRows, Real>
699 operator*(
700  Matrix<NumRows, NumCols, Real> const& M,
701  Vector<NumCols, Real> const& V)
702 {
703  Vector<NumRows, Real> result;
704  for (int r = 0; r < NumRows; ++r)
705  {
706  result[r] = (Real)0;
707  for (int c = 0; c < NumCols; ++c)
708  {
709  result[r] += M(r, c) * V[c];
710  }
711  }
712  return result;
713 }
714 
715 template <int NumRows, int NumCols, typename Real>
716 Vector<NumCols, Real> operator*(Vector<NumRows, Real> const& V,
717  Matrix<NumRows, NumCols, Real> const& M)
718 {
719  Vector<NumCols, Real> result;
720  for (int c = 0; c < NumCols; ++c)
721  {
722  result[c] = (Real)0;
723  for (int r = 0; r < NumRows; ++r)
724  {
725  result[c] += V[r] * M(r, c);
726  }
727  }
728  return result;
729 }
730 
731 template <int NumRows, int NumCols, int NumCommon, typename Real>
732 Matrix<NumRows, NumCols, Real>
733 operator*(
734  Matrix<NumRows, NumCommon, Real> const& A,
735  Matrix<NumCommon, NumCols, Real> const& B)
736 {
737  return MultiplyAB(A, B);
738 }
739 
740 template <int NumRows, int NumCols, int NumCommon, typename Real>
741 Matrix<NumRows, NumCols, Real>
742 MultiplyAB(
743  Matrix<NumRows, NumCommon, Real> const& A,
744  Matrix<NumCommon, NumCols, Real> const& B)
745 {
746  Matrix<NumRows, NumCols, Real> result;
747  for (int r = 0; r < NumRows; ++r)
748  {
749  for (int c = 0; c < NumCols; ++c)
750  {
751  result(r, c) = (Real)0;
752  for (int i = 0; i < NumCommon; ++i)
753  {
754  result(r, c) += A(r, i) * B(i, c);
755  }
756  }
757  }
758  return result;
759 }
760 
761 template <int NumRows, int NumCols, int NumCommon, typename Real>
762 Matrix<NumRows, NumCols, Real>
763 MultiplyABT(
764  Matrix<NumRows, NumCommon, Real> const& A,
765  Matrix<NumCols, NumCommon, Real> const& B)
766 {
767  Matrix<NumRows, NumCols, Real> result;
768  for (int r = 0; r < NumRows; ++r)
769  {
770  for (int c = 0; c < NumCols; ++c)
771  {
772  result(r, c) = (Real)0;
773  for (int i = 0; i < NumCommon; ++i)
774  {
775  result(r, c) += A(r, i) * B(c, i);
776  }
777  }
778  }
779  return result;
780 }
781 
782 template <int NumRows, int NumCols, int NumCommon, typename Real>
783 Matrix<NumRows, NumCols, Real>
784 MultiplyATB(
785  Matrix<NumCommon, NumRows, Real> const& A,
786  Matrix<NumCommon, NumCols, Real> const& B)
787 {
788  Matrix<NumRows, NumCols, Real> result;
789  for (int r = 0; r < NumRows; ++r)
790  {
791  for (int c = 0; c < NumCols; ++c)
792  {
793  result(r, c) = (Real)0;
794  for (int i = 0; i < NumCommon; ++i)
795  {
796  result(r, c) += A(i, r) * B(i, c);
797  }
798  }
799  }
800  return result;
801 }
802 
803 template <int NumRows, int NumCols, int NumCommon, typename Real>
804 Matrix<NumRows, NumCols, Real>
805 MultiplyATBT(
806  Matrix<NumCommon, NumRows, Real> const& A,
807  Matrix<NumCols, NumCommon, Real> const& B)
808 {
809  Matrix<NumRows, NumCols, Real> result;
810  for (int r = 0; r < NumRows; ++r)
811  {
812  for (int c = 0; c < NumCols; ++c)
813  {
814  result(r, c) = (Real)0;
815  for (int i = 0; i < NumCommon; ++i)
816  {
817  result(r, c) += A(i, r) * B(c, i);
818  }
819  }
820  }
821  return result;
822 }
823 
824 template <int NumRows, int NumCols, typename Real>
825 Matrix<NumRows, NumCols, Real>
826 MultiplyMD(
827  Matrix<NumRows, NumCols, Real> const& M,
828  Vector<NumCols, Real> const& D)
829 {
830  Matrix<NumRows, NumCols, Real> result;
831  for (int r = 0; r < NumRows; ++r)
832  {
833  for (int c = 0; c < NumCols; ++c)
834  {
835  result(r, c) = M(r, c) * D[c];
836  }
837  }
838  return result;
839 }
840 
841 template <int NumRows, int NumCols, typename Real>
842 Matrix<NumRows, NumCols, Real>
843 MultiplyDM(
844  Vector<NumRows, Real> const& D,
845  Matrix<NumRows, NumCols, Real> const& M)
846 {
847  Matrix<NumRows, NumCols, Real> result;
848  for (int r = 0; r < NumRows; ++r)
849  {
850  for (int c = 0; c < NumCols; ++c)
851  {
852  result(r, c) = D[r] * M(r, c);
853  }
854  }
855  return result;
856 }
857 
858 template <int NumRows, int NumCols, typename Real>
859 Matrix<NumRows, NumCols, Real>
860 OuterProduct(Vector<NumRows, Real> const& U, Vector<NumCols, Real> const& V)
861 {
862  Matrix<NumRows, NumCols, Real> result;
863  for (int r = 0; r < NumRows; ++r)
864  {
865  for (int c = 0; c < NumCols; ++c)
866  {
867  result(r, c) = U[r] * V[c];
868  }
869  }
870  return result;
871 }
872 
873 template <int N, typename Real>
874 void MakeDiagonal(Vector<N, Real> const& D, Matrix<N, N, Real>& M)
875 {
876  for (int i = 0; i < N*N; ++i)
877  {
878  M[i] = (Real)0;
879  }
880 
881  for (int i = 0; i < N; ++i)
882  {
883  M(i, i) = D[i];
884  }
885 }
886 
887 template <int N, typename Real>
888 Matrix<N + 1, N + 1, Real> HLift(Matrix<N, N, Real> const& M)
889 {
890  Matrix<N + 1, N + 1, Real> result;
891  result.MakeIdentity();
892  for (int r = 0; r < N; ++r)
893  {
894  for (int c = 0; c < N; ++c)
895  {
896  result(r, c) = M(r, c);
897  }
898  }
899  return result;
900 }
901 
902 // Extract the upper (N-1)-by-(N-1) block of the input N-by-N matrix.
903 template <int N, typename Real>
904 Matrix<N - 1, N - 1, Real> HProject(Matrix<N, N, Real> const& M)
905 {
906  static_assert(N >= 2, "Invalid matrix dimension.");
907  Matrix<N - 1, N - 1, Real> result;
908  for (int r = 0; r < N - 1; ++r)
909  {
910  for (int c = 0; c < N - 1; ++c)
911  {
912  result(r, c) = M(r, c);
913  }
914  }
915  return result;
916 }
917 
918 
919 // Matrix<N,C,Real>::Table
920 
921 template <int NumRows, int NumCols, typename Real> inline
922 Real const& Matrix<NumRows, NumCols, Real>::Table::operator()(int r, int c)
923  const
924 {
925 #if defined(GTE_USE_ROW_MAJOR)
926  return mStorage[r][c];
927 #else
928  return mStorage[c][r];
929 #endif
930 }
931 
932 template <int NumRows, int NumCols, typename Real> inline
933 Real& Matrix<NumRows, NumCols, Real>::Table::operator()(int r, int c)
934 {
935 #if defined(GTE_USE_ROW_MAJOR)
936  return mStorage[r][c];
937 #else
938  return mStorage[c][r];
939 #endif
940 }
941 
942 template <int NumRows, int NumCols, typename Real> inline
943 Real const& Matrix<NumRows, NumCols, Real>::Table::operator[](int i) const
944 {
945  Real const* elements = &mStorage[0][0];
946  return elements[i];
947 }
948 
949 template <int NumRows, int NumCols, typename Real> inline
950 Real& Matrix<NumRows, NumCols, Real>::Table::operator[](int i)
951 {
952  Real* elements = &mStorage[0][0];
953  return elements[i];
954 }
955 
956 
957 }
gte::Matrix::operator<=
bool operator<=(Matrix const &mat) const
Definition: GteMatrix.h:424
gte::Matrix::GetCol
Vector< NumRows, Real > GetCol(int c) const
Definition: GteMatrix.h:383
std::abs
gte::BSNumber< UIntegerType > abs(gte::BSNumber< UIntegerType > const &number)
Definition: GteBSNumber.h:966
gte::Matrix::operator!=
bool operator!=(Matrix const &mat) const
Definition: GteMatrix.h:412
gte::operator*=
DualQuaternion< Real > & operator*=(DualQuaternion< Real > &d, Real scalar)
gte::MakeDiagonal
void MakeDiagonal(GVector< Real > const &D, GMatrix< Real > &M)
Definition: GteGMatrix.h:829
gte::Matrix::MakeIdentity
void MakeIdentity()
Definition: GteMatrix.h:462
gte::MultiplyMD
GMatrix< Real > MultiplyMD(GMatrix< Real > const &M, GVector< Real > const &D)
Definition: GteGMatrix.h:781
gte::Matrix::Table::operator[]
Real const & operator[](int i) const
Definition: GteMatrix.h:943
gte::Matrix::operator>=
bool operator>=(Matrix const &mat) const
Definition: GteMatrix.h:436
values
GLenum GLsizei GLsizei GLint * values
Definition: glcorearb.h:1597
gte::MultiplyATB
GMatrix< Real > MultiplyATB(GMatrix< Real > const &A, GMatrix< Real > const &B)
Definition: GteGMatrix.h:737
gte::Matrix::operator==
bool operator==(Matrix const &mat) const
Definition: GteMatrix.h:406
gte::Matrix::operator[]
Real const & operator[](int i) const
Definition: GteMatrix.h:394
gte::Transpose
GMatrix< Real > Transpose(GMatrix< Real > const &M)
Definition: GteGMatrix.h:637
gte::HLift
GVector< Real > HLift(GVector< Real > const &v, Real last)
Definition: GteGVector.h:574
c
const GLubyte * c
Definition: glext.h:11671
gte::operator+
DualQuaternion< Real > operator+(DualQuaternion< Real > const &d)
gte::HProject
GVector< Real > HProject(GVector< Real > const &v)
Definition: GteGVector.h:587
gte::Matrix::GetRow
Vector< NumCols, Real > GetRow(int r) const
Definition: GteMatrix.h:372
gte::Matrix::Table::mStorage
std::array< std::array< Real, NumRows >, NumCols > mStorage
Definition: GteMatrix.h:103
gte::operator-=
DualQuaternion< Real > & operator-=(DualQuaternion< Real > &d0, DualQuaternion< Real > const &d1)
gte::Matrix::SetCol
void SetCol(int c, Vector< NumRows, Real > const &vec)
Definition: GteMatrix.h:362
gte::Matrix::operator>
bool operator>(Matrix const &mat) const
Definition: GteMatrix.h:430
gte::MultiplyAB
GMatrix< Real > MultiplyAB(GMatrix< Real > const &A, GMatrix< Real > const &B)
Definition: GteGMatrix.h:693
gte::Matrix::operator<
bool operator<(Matrix const &mat) const
Definition: GteMatrix.h:418
gte::Matrix::operator()
Real const & operator()(int r, int c) const
Definition: GteMatrix.h:340
r
GLboolean r
Definition: glcorearb.h:1217
gte::Matrix::Unit
static Matrix Unit(int r, int c)
Definition: GteMatrix.h:481
gte::MultiplyATBT
GMatrix< Real > MultiplyATBT(GMatrix< Real > const &A, GMatrix< Real > const &B)
Definition: GteGMatrix.h:759
gte::LInfinityNorm
Real LInfinityNorm(GMatrix< Real > const &M)
Definition: GteGMatrix.h:575
gte::Matrix::MakeZero
void MakeZero()
Definition: GteMatrix.h:442
gte::Matrix::SetRow
void SetRow(int r, Vector< NumCols, Real > const &vec)
Definition: GteMatrix.h:352
gte::MultiplyABT
GMatrix< Real > MultiplyABT(GMatrix< Real > const &A, GMatrix< Real > const &B)
Definition: GteGMatrix.h:715
gte::operator-
DualQuaternion< Real > operator-(DualQuaternion< Real > const &d)
gte::operator+=
DualQuaternion< Real > & operator+=(DualQuaternion< Real > &d0, DualQuaternion< Real > const &d1)
gte::Inverse
Quaternion< Real > Inverse(Quaternion< Real > const &d)
Definition: GteQuaternion.h:227
gte::Matrix::Identity
static Matrix Identity()
Definition: GteMatrix.h:490
gte::L2Norm
Real L2Norm(GMatrix< Real > const &M)
Definition: GteGMatrix.h:564
gte::Determinant
Real Determinant(GMatrix< Real > const &M)
Definition: GteGMatrix.h:618
gte::operator/=
DualQuaternion< Real > & operator/=(DualQuaternion< Real > &d, Real scalar)
gte::operator*
Vector4< float > operator*(Transform const &M, Vector4< float > const &V)
Definition: GteTransform.cpp:487
result
GLuint64EXT * result
Definition: glext.h:10003
gte::Matrix::Table::operator()
Real const & operator()(int r, int c) const
Definition: GteMatrix.h:922
gte::Matrix::MakeUnit
void MakeUnit(int r, int c)
Definition: GteMatrix.h:452
gte::OuterProduct
GMatrix< Real > OuterProduct(GVector< Real > const &U, GVector< Real > const &V)
Definition: GteGMatrix.h:815
gte::L1Norm
Real L1Norm(GMatrix< Real > const &M)
Definition: GteGMatrix.h:553
gte::operator/
DualQuaternion< Real > operator/(DualQuaternion< Real > const &d, Real scalar)
gte::MultiplyDM
GMatrix< Real > MultiplyDM(GVector< Real > const &D, GMatrix< Real > const &M)
Definition: GteGMatrix.h:798
gte::Matrix::Zero
static Matrix Zero()
Definition: GteMatrix.h:473
gte::Matrix::mTable
Table mTable
Definition: GteMatrix.h:107

geometric_tools_engine
Author(s): Yijiang Huang
autogenerated on Thu Jul 18 2019 04:00:01