Simbody  3.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MatrixCharacteristics.h
Go to the documentation of this file.
1 #ifndef SimTK_SIMMATRIX_MATRIX_CHARACTERISTICS_H_
2 #define SimTK_SIMMATRIX_MATRIX_CHARACTERISTICS_H_
3 
4 /* -------------------------------------------------------------------------- *
5  * Simbody(tm): SimTKcommon *
6  * -------------------------------------------------------------------------- *
7  * This is part of the SimTK biosimulation toolkit originating from *
8  * Simbios, the NIH National Center for Physics-Based Simulation of *
9  * Biological Structures at Stanford, funded under the NIH Roadmap for *
10  * Medical Research, grant U54 GM072970. See https://simtk.org/home/simbody. *
11  * *
12  * Portions copyright (c) 2005-13 Stanford University and the Authors. *
13  * Authors: Michael Sherman *
14  * Contributors: *
15  * *
16  * Licensed under the Apache License, Version 2.0 (the "License"); you may *
17  * not use this file except in compliance with the License. You may obtain a *
18  * copy of the License at http://www.apache.org/licenses/LICENSE-2.0. *
19  * *
20  * Unless required by applicable law or agreed to in writing, software *
21  * distributed under the License is distributed on an "AS IS" BASIS, *
22  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. *
23  * See the License for the specific language governing permissions and *
24  * limitations under the License. *
25  * -------------------------------------------------------------------------- */
26 
33 #include "SimTKcommon/Scalar.h"
34 
35 #include <iostream>
36 #include <cassert>
37 #include <complex>
38 #include <cstddef>
39 #include <utility> // for std::pair
40 
41 namespace SimTK {
42 
43 
44 class MatrixStructure;
45 class MatrixStorage;
46 class MatrixOutline;
47 class MatrixCondition;
48 class MatrixCharacter;
49 class MatrixCommitment;
50 
51 
52 // ------------------------------ MatrixStructure -----------------------------
57 // ----------------------------------------------------------------------------
58 class SimTK_SimTKCOMMON_EXPORT MatrixStructure {
59 public:
60  enum Structure {
61  NoStructure = 0x00000000,
62  Matrix1d = 0x00000001,
63  Zero = 0x00000002,
64  Identity = 0x00000004,
65  Permutation = 0x00000008,
66  RepeatedDiagonal = 0x00000010,
67  Diagonal = 0x00000020,
68  BiDiagonal = 0x00000040,
69  TriDiagonal = 0x00000080,
70  BandedSymmetric = 0x00000100,
71  BandedHermitian = 0x00000200,
72  Banded = 0x00000400,
73  Triangular = 0x00000800,
74  QuasiTriangular = 0x00001000,
75  Hessenberg = 0x00002000,
76  Symmetric = 0x00004000,
77  Hermitian = 0x00008000,
78  SkewSymmetric = 0x00010000,
79  SkewHermitian = 0x00020000,
80  Full = 0x00040000
81  };
82  static const char* name(Structure);
83 
84  typedef unsigned int StructureMask; // 32 bits
85  static const StructureMask AnyStructure = 0x0007ffffU; // see above
86  static const StructureMask UncommittedStructure = 0xffffffffU;
87  static StructureMask calcStructureMask(Structure);
88 
93  enum Position {
94  NoPosition = 0x0000,
95  Lower = 0x0001, // matrix is lower triangular (default)
96  Upper = 0x0002 // matrix is upper triangular
97  };
98  static const char* name(Position);
99 
100  typedef unsigned short PositionMask; // 16 bits
101  static const PositionMask AnyPosition = 0x0003U; // see above
102  static const PositionMask UncommittedPosition = 0xffffU;
103  static PositionMask calcPositionMask(Structure);
104 
110  enum DiagValue {
111  NoDiagValue = 0x0000,
112  StoredDiag = 0x0001, // could be anything (default)
113  ZeroDiag = 0x0002, // zero (e.g. for skew matrices)
114  UnitDiag = 0x0004 // unit (one) diagonal is used frequently by Lapack
115  };
116  static const char* name(DiagValue);
117 
118  typedef unsigned short DiagValueMask; // 16 bits
119  static const DiagValueMask AnyDiagValue = 0x0003U;
120  static const DiagValueMask UncommittedDiagValue = 0xffffU;
121  static DiagValueMask calcDiagValueMask(Structure);
122 
123  MatrixStructure& setMissingAttributes() {
124  if (structure == NoStructure)
125  structure = Full;
126  if (position == NoPosition)
127  position = Lower;
128  if (diagValue == NoDiagValue)
129  diagValue = StoredDiag;
130  return *this;
131  }
132 
133  std::string name() const {
134  return std::string(name(getStructure()))
135  + "|" + std::string(name(getPosition()))
136  + "|" + std::string(name(getDiagValue()));
137  }
138 
139  struct Mask {
140  Mask() {setToUncommitted();}
141  Mask(StructureMask sm, PositionMask pm, DiagValueMask dm)
142  : structure(sm), position(pm), diagValue(dm) {}
143  Mask& setToUncommitted()
144  { structure=UncommittedStructure; position=UncommittedPosition;
145  diagValue=UncommittedDiagValue; return *this; }
146  bool isUncommitted() const
147  { return structure==UncommittedStructure && position==UncommittedPosition
148  && diagValue==UncommittedDiagValue; }
149  bool isSatisfiedBy(Structure str, Position pos, DiagValue diag) const
150  { return ((StructureMask)str&structure)==(StructureMask)str
151  && ((PositionMask)pos&position)==(PositionMask)pos
152  && ((DiagValueMask)diag&diagValue)==(DiagValueMask)diag; }
153  bool isSatisfiedBy(const MatrixStructure& actual) const
154  { return isSatisfiedBy(actual.getStructure(), actual.getPosition(),
155  actual.getDiagValue()); }
156 
157  StructureMask structure;
158  PositionMask position;
159  DiagValueMask diagValue;
160  };
161 
162  MatrixStructure() {setToNone();}
163 
166  MatrixStructure(Structure s, Position p=NoPosition, DiagValue d=NoDiagValue)
167  : structure(s), position(p), diagValue(d) {}
168 
173  Mask mask() const;
174 
175  Structure getStructure() const {return structure;}
176  Position getPosition() const {return position;}
177  DiagValue getDiagValue() const {return diagValue;}
178 
179  MatrixStructure& setStructure(Structure s) {structure=s; return *this;}
180  MatrixStructure& setPosition (Position p) {position=p; return *this;}
181  MatrixStructure& setDiagValue(DiagValue d) {diagValue=d; return *this;}
182 
183  MatrixStructure& set(Structure s, Position p, DiagValue d)
184  { structure=s; position=p; diagValue=d; return *this; }
185 
186  MatrixStructure& setToNone()
187  { structure=NoStructure; position=NoPosition;
188  diagValue=NoDiagValue; return *this; }
189 
190 private:
191  Structure structure:32;
192  Position position:16;
193  DiagValue diagValue:16;
194 };
195 
196 
197 // ------------------------------ MatrixStorage -------------------------------
202 // ----------------------------------------------------------------------------
203 class SimTK_SimTKCOMMON_EXPORT MatrixStorage {
204 public:
205  enum Packing {
206  NoPacking = 0x0000,
207  Full = 0x0001, // full storage layout
208  TriInFull = 0x0002, // a triangular piece of a full storage layout
209  TriPacked = 0x0004, // triangle packed into minimal storage, at performance cost
210  Banded = 0x0008, // a packed, banded storage format
211  Vector = 0x0010, // a possibly-strided or scattered vector
212  Scalar = 0x0020, // a single scalar is stored
213  Permutation = 0x0040 // a permuted identity matrix
214  };
215  static const char* name(Packing);
216  typedef unsigned short PackingMask;
217  static const PackingMask AllPacking = 0x007fU; // see above
218  static const PackingMask UncommittedPacking = 0xffffU;
219 
220  enum Placement {
221  NoPlacement = 0x0000,
222  Lower = 0x0001, // stored in lower triangle of full matrix
223  Upper = 0x0002, // stored in upper triangle of full matrix
224  };
225  static const char* name(Placement);
226  typedef unsigned short PlacementMask;
227  static const PlacementMask AllPlacement = 0x0003U; // see above
228  static const PlacementMask UncommittedPlacement = 0xffffU;
229 
230  enum Order {
231  NoOrder = 0x0000,
232  ColumnOrder = 0x0001, // matrix is stored by columns
233  RowOrder = 0x0002, // matrix is stored by rows
234  };
235  static const char* name(Order);
236  typedef unsigned short OrderMask;
237  static const OrderMask AllOrder = 0x03U; // see above
238  static const OrderMask UncommittedOrder = 0xffU;
239 
240  enum Diagonal {
241  NoDiag = 0x0000,
242  StoredDiag = 0x0001, // matrix diagonal is stored
243  AssumedDiag = 0x0002 // matrix diagonal is not stored but has known value
244  };
245  static const char* name(Diagonal);
246  typedef unsigned short DiagonalMask;
247  static const DiagonalMask AllDiagonal = 0x0003U; // see above
248  static const DiagonalMask UncommittedDiagonal = 0xffffU;
249 
252  struct Mask {
253  Mask()
254  : packing(UncommittedPacking), placement(UncommittedPlacement),
255  order(UncommittedOrder), diagonal(UncommittedDiagonal) {}
256  Mask(PackingMask pkm, PlacementMask plm, OrderMask om, DiagonalMask dm)
257  : packing(pkm), placement(plm), order(om), diagonal(dm) {}
258  Mask& setToUncommitted()
259  { packing=UncommittedPacking; placement=UncommittedPlacement;
260  order=UncommittedOrder; diagonal=UncommittedDiagonal; return *this; }
261  bool isUncommitted() const
262  { return packing==UncommittedPacking && placement==UncommittedPlacement
263  && order==UncommittedOrder && diagonal==UncommittedDiagonal; }
264  bool isSatisfiedBy(Packing pack, Placement place, Order ord, Diagonal diag) const
265  { return ((PackingMask)pack & packing) == (PackingMask) pack
266  && ((PlacementMask)place & placement) == (PlacementMask)place
267  && ((OrderMask)ord & order) == (OrderMask) ord
268  && ((DiagonalMask)diag & diagonal) == (DiagonalMask) diag; }
269  bool isSatisfiedBy(const MatrixStorage& actual) const
270  { return isSatisfiedBy(actual.getPacking(), actual.getPlacement(),
271  actual.getOrder(), actual.getDiagonal());}
272 
273  PackingMask packing;
274  PlacementMask placement;
275  OrderMask order;
276  DiagonalMask diagonal;
277  };
278 
279  static MatrixStorage calcDefaultStorage(const MatrixStructure&,
280  const MatrixOutline&);
281 
282  std::string name() const {
283  return std::string(name(getPacking()))
284  + "|" + std::string(name(getPlacement()))
285  + "|" + std::string(name(getOrder()))
286  + "|" + std::string(name(getDiagonal()));
287  }
288 
293  Mask mask() const {
294  Mask ms; // initially uncommitted
295  if (packing) ms.packing = (PackingMask)packing;
296  if (placement) ms.placement = (PlacementMask)placement;
297  if (order) ms.order = (OrderMask)order;
298  if (diagonal) ms.diagonal = (DiagonalMask)diagonal;
299  return ms;
300  }
301 
303  MatrixStorage()
304  : packing(NoPacking), placement(NoPlacement), order(NoOrder), diagonal(NoDiag) {}
305 
309  MatrixStorage(Packing pk, Placement pl=NoPlacement, Order o=NoOrder, Diagonal d=NoDiag)
310  : packing(pk), placement(pl), order(o), diagonal(d) {}
311 
314  MatrixStorage(Packing pk, Order o)
315  : packing(pk), placement(NoPlacement), order(o), diagonal(StoredDiag) {}
316 
319  MatrixStorage& setMissingAttributes() {
320  if (packing==NoPacking)
321  packing = Full;
322  if (placement==NoPlacement)
323  placement = Lower;
324  if (order==NoOrder)
325  order = ColumnOrder;
326  if (diagonal==NoDiag)
327  diagonal = StoredDiag;
328  return *this;
329  }
330 
332  MatrixStorage& setToNone()
333  { packing=NoPacking; placement=NoPlacement;
334  order=NoOrder; diagonal=NoDiag; return *this; }
335 
336  MatrixStorage& setPacking(Packing p) {packing = p; return *this;}
337  MatrixStorage& setPlacement(Placement p) {placement = p; return *this;}
338  MatrixStorage& setOrder(Order o) {order = o; return *this;}
339  MatrixStorage& setDiagonal(Diagonal d) {diagonal = d; return *this;}
340 
341  Packing getPacking() const {return packing;}
342  Placement getPlacement() const {return placement;}
343  Order getOrder() const {return order;}
344  Diagonal getDiagonal() const {return diagonal;}
345 
346 private:
347  Packing packing:16;
348  Placement placement:16;
349  Order order:16;
350  Diagonal diagonal:16;
351 };
352 
353 
354 // ------------------------------- MatrixOutline ------------------------------
375 // ----------------------------------------------------------------------------
376 class SimTK_SimTKCOMMON_EXPORT MatrixOutline {
377 public:
378  enum Outline {
379  NoOutline = 0x0000,
380  Scalar = 0x0001, // 1x1
381  Column = 0x0002, // mx1, m != 1
382  Row = 0x0004, // 1xn, n != 1
383  Square = 0x0008, // mxn, m == n
384  Wide = 0x0010, // mxn, m < n
385  Tall = 0x0020, // mxn, m > n
386  Rectangular = 0x0040 // mxn
387  };
388  static const char* name(Outline);
389 
390  typedef unsigned short OutlineMask;
391  static const OutlineMask AnyOutline = 0x007fU; // see above
392  static const OutlineMask UncommittedOutline = 0xffffU;
393 
394  struct Mask {
395  Mask() : outline(UncommittedOutline) {}
396  explicit Mask(OutlineMask mask) : outline(mask) {}
397  Mask& setToUncommitted() {outline=UncommittedOutline; return *this;}
398  bool isUncommitted() const {return outline==UncommittedOutline;}
399  bool isSatisfiedBy(const MatrixOutline& actual) const
400  { return ((OutlineMask)actual.outline & outline) == (OutlineMask)actual.outline; }
401 
402  OutlineMask outline;
403  };
404 
405  std::string name() const {return std::string(name(getOutline()));}
406 
409  MatrixOutline() : outline(NoOutline) {}
410 
412  MatrixOutline(Outline outline) : outline(outline) {}
413 
415  MatrixOutline& setToNone() {outline=NoOutline; return *this;}
416 
420  static OutlineMask calcMask(Outline);
421 
424  Mask mask() const {return Mask(calcMask(getOutline()));}
425 
427  bool isSizeOK(int m, int n) const;
428 
430  void getMinimumSize(int& m, int& n) const;
431 
433  static MatrixOutline calcFromSize(int m, int n);
434 
436  Outline getOutline() const {return outline;}
437 
438 private:
439  Outline outline:16;
440 };
441 
442 
443 
444 // ---------------------------- MatrixCondition -------------------------------
453 // ----------------------------------------------------------------------------
454 class SimTK_SimTKCOMMON_EXPORT MatrixCondition {
455 public:
456  enum Condition {
457  UnknownCondition = 0x0000,
458  Orthogonal = 0x0001, // implies well conditioned
459  PositiveDefinite = 0x0002, // implies well conditioned
460  WellConditioned = 0x0004, // implies full rank
461  FullRank = 0x0008, // but might have bad conditioning
462  Singular = 0x0010 // implies possible bad conditioning
463  };
464  static const char* name(Condition);
465 
466  typedef unsigned short ConditionMask; // 16 bits in mask
467  static const ConditionMask AnyCondition = 0x001fU; // see above
468  static const ConditionMask UncommittedCondition = 0xffffU;
469 
470  enum Diagonal {
471  UnknownDiagonal = 0x0000,
472  ZeroDiagonal = 0x0001,
473  OneDiagonal = 0x0002,
474  RealDiagonal = 0x0004,
475  ImaginaryDiagonal = 0x0008
476  };
477  static const char* name(Diagonal);
478 
479  typedef unsigned short DiagonalMask; // 16 bits in mask
480  static const DiagonalMask AnyDiagonal = 0x000fU; // see above
481  static const DiagonalMask UncommittedDiagonal = 0xffffU;
482 
484  struct Mask {
485  Mask() : condition(UncommittedCondition), diagonal(UncommittedDiagonal) {}
486  Mask(ConditionMask cmask, DiagonalMask dmask) : condition(cmask), diagonal(dmask) {}
487  Mask& setToUncommitted()
488  { condition=UncommittedCondition; diagonal=UncommittedDiagonal; return *this;}
489  bool isUncommitted() const
490  { return condition==UncommittedCondition && diagonal==UncommittedDiagonal;}
491  bool isSatisfiedBy(const MatrixCondition& actual) const
492  { return ((ConditionMask)actual.condition & condition) == (ConditionMask)actual.condition
493  && ((DiagonalMask) actual.diagonal & diagonal) == (DiagonalMask)actual.diagonal; }
494 
495  ConditionMask condition;
496  DiagonalMask diagonal;
497  };
498 
499  std::string name() const
500  { return std::string(name(getCondition())) + "|" + std::string(name(getDiagonal()));}
501 
504  MatrixCondition() : condition(UnknownCondition), diagonal(UnknownDiagonal) {}
505 
508  MatrixCondition(Condition cond, Diagonal diag=UnknownDiagonal) : condition(cond) {}
509 
511  MatrixCondition& setToNone() {condition=UnknownCondition; diagonal=UnknownDiagonal; return *this;}
512 
518  static ConditionMask calcMask(Condition);
519 
525  static DiagonalMask calcMask(Diagonal);
526 
529  Mask mask() const
530  { return Mask(calcMask(getCondition()), calcMask(getDiagonal())); }
531 
532  Condition getCondition() const {return condition;}
533  Diagonal getDiagonal() const {return diagonal;}
534 
535  MatrixCondition& setCondition(Condition c) {condition=c; return *this;}
536  MatrixCondition& setDiagonal (Diagonal d) {diagonal=d; return *this;}
537 
538 private:
539  Condition condition:16;
540  Diagonal diagonal:16;
541 };
542 
543 
544 
545 // ------------------------------ MatrixCharacter -----------------------------
596 class SimTK_SimTKCOMMON_EXPORT MatrixCharacter {
597 public:
600  MatrixCharacter() : nr(0), nc(0), lband(0), uband(0) {}
601 
602  // Some handy predefined MatrixCharacters.
603  class LapackFull;
604  class Vector;
605  class RowVector;
606 
608  MatrixCharacter& setToNone() {
609  nr=nc=lband=uband=0;
610  structure.setToNone(); outline.setToNone();
611  storage.setToNone(); condition.setToNone();
612  return *this;
613  }
614 
617  int nrow() const {return nr;}
618  int ncol() const {return nc;}
619  std::pair<int,int> getSize() const {return std::pair<int,int>(nrow(),ncol());}
620  ptrdiff_t nelt() const {return (ptrdiff_t)nrow() * (ptrdiff_t)ncol();}
621 
622  int getLowerBandwidth() const {return lband;}
623  int getUpperBandwidth() const {return uband;}
624  std::pair<int,int> getBandwidth() const
625  { return std::pair<int,int>(getLowerBandwidth(), getUpperBandwidth()); }
626 
627  const MatrixStructure& getStructure() const {return structure;}
628  const MatrixStorage& getStorage() const {return storage;}
629  const MatrixOutline& getOutline() const {return outline;}
630  const MatrixCondition& getCondition() const {return condition;}
631 
632  MatrixStructure& updStructure() {return structure;}
633  MatrixStorage& updStorage() {return storage;}
634  MatrixOutline& updOutline() {return outline;}
635  MatrixCondition& updCondition() {return condition;}
636 
637  MatrixCharacter& setStructure(const MatrixStructure& sa) {structure = sa; return *this;}
638  MatrixCharacter& setStorage (const MatrixStorage& sa) {storage = sa; return *this;}
639  MatrixCharacter& setOutline (const MatrixOutline& oa) {outline = oa; return *this;}
640  MatrixCharacter& setCondition(const MatrixCondition& ca) {condition = ca; return *this;}
641 
642 
644  MatrixCharacter& setActualSize(int m, int n)
645  { setSize(m,n); outline = MatrixOutline::calcFromSize(m,n); return *this; }
646  MatrixCharacter& setActualNumRows(int m)
647  { setNumRows(m); outline = MatrixOutline::calcFromSize(m,ncol()); return *this; }
648  MatrixCharacter& setActualNumCols(int n)
649  { setNumCols(n); outline = MatrixOutline::calcFromSize(nrow(),n); return *this; }
650 
651  MatrixCharacter& setBandwidth(int lb, int ub) {
652  assert(lb>=0 && lb>=0);
653  lband = lb; uband = ub;
654  return *this;
655  }
656  MatrixCharacter& setLowerBandwidth(int lb) {
657  assert(lb>=0);
658  lband = lb;
659  return *this;
660  }
661  MatrixCharacter& setUpperBandwidth(int ub) {
662  assert(ub>=0);
663  uband = ub;
664  return *this;
665  }
666 
667  class Mask; // defined below
668 
669 protected:
670  MatrixCharacter(int m, int n,
671  int lb, int ub,
672  MatrixStructure structure,
673  MatrixStorage storage,
674  MatrixCondition condition)
675  : nr(m), nc(n), lband(lb), uband(ub),
676  structure(structure), storage(storage),
677  outline(MatrixOutline::calcFromSize(m,n)),
678  condition(condition) {}
679 
680 
681  int nr,
682  nc;
683  int lband,
684  uband;
685  MatrixStructure structure;
686  MatrixStorage storage;
687  MatrixOutline outline;
688  MatrixCondition condition;
689 
690 private:
691  // These are private because they don't set the outline as well.
692  MatrixCharacter& setSize(int m, int n)
693  { assert(m>=0 && n>=0); nr = m; nc = n; return *this; }
694  MatrixCharacter& setNumRows(int m)
695  { assert(m>=0); nr = m; return *this; }
696  MatrixCharacter& setNumCols(int n)
697  { assert(n>=0); nc = n; return *this; }
698 };
699 
702 SimTK_SimTKCOMMON_EXPORT std::ostream&
703 operator<<(std::ostream& o, const MatrixCharacter&);
704 
710 class MatrixCharacter::LapackFull : public MatrixCharacter {
711 public:
712  LapackFull(int m, int n)
713  : MatrixCharacter(m,n,0,0,
714  MatrixStructure(MatrixStructure::Full),
715  MatrixStorage(MatrixStorage::Full,MatrixStorage::ColumnOrder),
716  MatrixCondition()) {}
717 };
718 
723 class MatrixCharacter::Vector : public MatrixCharacter {
724 public:
725  Vector(int m)
726  : MatrixCharacter(m,1,0,0,
727  MatrixStructure(MatrixStructure::Matrix1d),
728  MatrixStorage(MatrixStorage::Vector,MatrixStorage::ColumnOrder),
729  MatrixCondition()) {}
730 };
731 
736 class MatrixCharacter::RowVector : public MatrixCharacter {
737 public:
738  RowVector(int n)
739  : MatrixCharacter(1,n,0,0,
740  MatrixStructure(MatrixStructure::Matrix1d),
741  MatrixStorage(MatrixStorage::Vector,MatrixStorage::RowOrder),
742  MatrixCondition()) {}
743 };
744 
745 // -------------------------- MatrixCharacter::Mask ---------------------------
748 // ----------------------------------------------------------------------------
749 class MatrixCharacter::Mask {
750 public:
751  Mask() {setToUncommitted();}
752 
753  typedef unsigned int SizeMask;
754  static const SizeMask SizeUncommitted = 0xffffffffU;
755 
756  bool isResizeable() const {return nr==SizeUncommitted || nc==SizeUncommitted;}
757  bool isFullyResizeable() const {return nr==SizeUncommitted && nc==SizeUncommitted;}
758  bool isNumRowsLocked() const {return nr!=SizeUncommitted;}
759  bool isNumColsLocked() const {return nc!=SizeUncommitted;}
760 
761  unsigned int getNumRowsMask() const {return nr;}
762  unsigned int getNumColsMask() const {return nc;}
763  unsigned int getLowerBandwidthMask() const {return lband;}
764  unsigned int getUpperBandwidthMask() const {return uband;}
765 
766  int getDefaultNumRows() const {return isNumRowsLocked() ? nr : 0;}
767  int getDefaultNumCols() const {return isNumColsLocked() ? nc : 0;}
768 
769  bool isLowerBandwidthLocked() const {return lband!=SizeUncommitted;}
770  bool isUpperBandwidthLocked() const {return uband!=SizeUncommitted;}
771  int getDefaultLowerBandwidth() const {return isLowerBandwidthLocked() ? lband : 0;}
772  int getDefaultUpperBandwidth() const {return isUpperBandwidthLocked() ? uband : 0;}
773 
775  Mask& setToUncommitted() {
776  nr=nc=lband=uband=SizeUncommitted;
777  structure.setToUncommitted(); storage.setToUncommitted();
778  outline.setToUncommitted(); condition.setToUncommitted();
779  return *this;
780  }
781 
783  bool isUncommitted() const {
784  return nr==SizeUncommitted && nc==SizeUncommitted
785  && lband==SizeUncommitted && uband==SizeUncommitted
786  && structure.isUncommitted() && storage.isUncommitted()
787  && outline.isUncommitted() && condition.isUncommitted();
788  }
789 
791  bool isSatisfiedBy(const MatrixCharacter& actual) const {
792  return isSizeOK(actual.nr, actual.nc)
793  && isBandwidthOK(actual.lband, actual.uband)
794  && structure.isSatisfiedBy(actual.getStructure())
795  && storage.isSatisfiedBy(actual.getStorage())
796  && outline.isSatisfiedBy(actual.getOutline())
797  && condition.isSatisfiedBy(actual.getCondition());
798  }
799 
801  bool isSizeOK(int m, int n) const
802  { return ((SizeMask)m & nr) == (SizeMask)m
803  && ((SizeMask)n & nc) == (SizeMask)n; }
804 
807  bool isBandwidthOK(int lower, int upper) const
808  { return ((SizeMask)lower & lband) == (SizeMask)lower
809  && ((SizeMask)upper & uband) == (SizeMask)upper; }
810 
811  SizeMask nr,
812  nc;
813  SizeMask lband,
814  uband;
815  MatrixStructure::Mask structure;
816  MatrixStorage::Mask storage;
817  MatrixOutline::Mask outline;
818  MatrixCondition::Mask condition;
819 
820 friend class MatrixCommitment;
821 };
822 
823 // ----------------------------- MatrixCommitment -----------------------------
824 
829 
830 // ----------------------------------------------------------------------------
831 class SimTK_SimTKCOMMON_EXPORT MatrixCommitment {
832 public:
833  MatrixCommitment() {} // set commitments to "none" and masks to "uncommitted"
834 
837  MatrixCommitment(const MatrixStructure& str)
838  { new (this) MatrixCommitment(str, MatrixStorage(), MatrixOutline(), MatrixCondition());}
839 
840  class Vector;
841  class RowVector;
842  class Triangular;
843  class Symmetric;
844  class Hermitian;
845  class SkewSymmetric;
846  class SkewHermitian;
847 
848  MatrixCommitment& commitSize(int m, int n)
849  { commitNumRows(m); commitNumCols(n); return *this; }
850  MatrixCommitment& commitNumRows(int m)
851  { SimTK_SIZECHECK_NONNEG(m, "MatrixCommitment::commitNumRows()");
852  masks.nr = m; return *this; }
853  MatrixCommitment& commitNumCols(int n)
854  { SimTK_SIZECHECK_NONNEG(n, "MatrixCommitment::commitNumCols()");
855  masks.nc = n; return *this; }
856 
857  MatrixCommitment& commitBandwidth(int lb, int ub)
858  { commitLowerBandwidth(lb); commitUpperBandwidth(ub); return *this;}
859  MatrixCommitment& commitLowerBandwidth(int lb)
860  { SimTK_SIZECHECK_NONNEG(lb, "MatrixCommitment::commitLowerBandwidth()");
861  masks.lband = lb; return *this; }
862  MatrixCommitment& commitUpperBandwidth(int ub)
863  { SimTK_SIZECHECK_NONNEG(ub, "MatrixCommitment::commitUpperBandwidth()");
864  masks.uband = ub; return *this; }
865 
866  MatrixCommitment& commitStructure(const MatrixStructure& s)
867  { structure=s; masks.structure=s.mask(); return *this; }
868  MatrixCommitment& commitStorage (const MatrixStorage& s)
869  { storage=s; masks.storage =s.mask(); return *this; }
870  MatrixCommitment& commitOutline (const MatrixOutline& o)
871  { outline=o; masks.outline =o.mask(); return *this; }
872  MatrixCommitment& commitCondition(const MatrixCondition& c)
873  { condition=c; masks.condition=c.mask(); return *this; }
874 
885  MatrixCharacter calcDefaultCharacter(int minNumRows, int minNumCols) const;
886 
888  const MatrixStructure& getStructureCommitment() const {return structure;}
889  const MatrixStorage& getStorageCommitment() const {return storage;}
890  const MatrixOutline& getOutlineCommitment() const {return outline;}
891  const MatrixCondition& getConditionCommitment() const {return condition;}
892 
894  const MatrixStructure::Mask& getStructureMask() const {return masks.structure;}
895  const MatrixStorage::Mask& getStorageMask() const {return masks.storage;}
896  const MatrixOutline::Mask& getOutlineMask() const {return masks.outline;}
897  const MatrixCondition::Mask& getConditionMask() const {return masks.condition;}
898 
899  MatrixCharacter::Mask::SizeMask getNumRowsMask() const {return masks.nr;}
900  MatrixCharacter::Mask::SizeMask getNumColsMask() const {return masks.nc;}
901  MatrixCharacter::Mask::SizeMask getLowerBandwidthMask() const {return masks.lband;}
902  MatrixCharacter::Mask::SizeMask getUpperBandwidthMask() const {return masks.uband;}
903 
904  int getDefaultNumRows() const {return masks.getDefaultNumRows();}
905  int getDefaultNumCols() const {return masks.getDefaultNumRows();}
906 
907  bool isSizeOK(int m, int n) const {return masks.isSizeOK(m,n);}
908  bool isSizeOK(const std::pair<int,int>& mn) const
909  { return isSizeOK(mn.first, mn.second); }
910 
911  bool isBandwidthOK(int lower, int upper) const {return masks.isBandwidthOK(lower,upper);}
912 
913  bool isSatisfiedBy(const MatrixCharacter& actual) const
914  { return masks.isSatisfiedBy(actual); }
915  bool isStructureOK(const MatrixStructure& s) const
916  { return getStructureMask().isSatisfiedBy(s); }
917  bool isStorageOK(const MatrixStorage& s) const
918  { return getStorageMask().isSatisfiedBy(s); }
919  bool isOutlineOK(const MatrixOutline& o) const
920  { return getOutlineMask().isSatisfiedBy(o); }
921  bool isConditionOK(const MatrixCondition& c) const
922  { return getConditionMask().isSatisfiedBy(c); }
923 
924  bool isResizeable() const {return masks.isResizeable();}
925  bool isFullyResizeable() const {return masks.isFullyResizeable();;}
926  bool isNumRowsLocked() const {return masks.isNumRowsLocked();}
927  bool isNumColsLocked() const {return masks.isNumColsLocked();}
928 
929  bool isStructureCommitted() const
930  { return !getStructureMask().isUncommitted(); }
931  bool isStorageCommitted() const
932  { return !getStorageMask().isUncommitted();}
933  bool isOutlineCommitted() const
934  { return !getOutlineMask().isUncommitted(); }
935  bool isConditionCommitted() const
936  { return !getConditionMask().isUncommitted();}
937 
939  void clear() {
940  structure.setToNone();
941  storage.setToNone();
942  outline.setToNone();
943  condition.setToNone();
944  masks.setToUncommitted();
945  }
946 
947 protected:
948  MatrixCommitment(const MatrixStructure& structure,
949  const MatrixStorage& storage,
950  const MatrixOutline& outline,
951  const MatrixCondition& condition)
952  : structure(structure), storage(storage),
953  outline(outline), condition(condition),
954  masks() // set to all 1's
955  {
956  if (outline.getOutline()==MatrixOutline::Scalar) commitSize(1,1);
957  else if (outline.getOutline()==MatrixOutline::Column) commitNumCols(1);
958  else if (outline.getOutline()==MatrixOutline::Row) commitNumRows(1);
959 
960  masks.structure = structure.mask();
961  masks.storage = storage.mask();
962  masks.outline = outline.mask();
963  masks.condition = condition.mask();
964  }
965 
968  MatrixStructure structure;
969  MatrixStorage storage;
970  MatrixOutline outline;
971  MatrixCondition condition;
972 
975  MatrixCharacter::Mask masks;
976 };
977 
978 
980 class MatrixCommitment::Vector : public MatrixCommitment {
981 public:
983  Vector()
984  : MatrixCommitment
985  ( MatrixStructure(MatrixStructure::Matrix1d),
986  MatrixStorage(),
987  MatrixOutline(MatrixOutline::Column),
988  MatrixCondition())
989  {
990  }
992  explicit Vector(int m)
993  : MatrixCommitment
994  ( MatrixStructure(MatrixStructure::Matrix1d),
995  MatrixStorage(),
996  MatrixOutline(MatrixOutline::Column),
997  MatrixCondition())
998  {
999  commitNumRows(m);
1000  }
1001 };
1002 
1004 class MatrixCommitment::RowVector : public MatrixCommitment {
1005 public:
1007  RowVector()
1008  : MatrixCommitment
1009  ( MatrixStructure(MatrixStructure::Matrix1d),
1010  MatrixStorage(),
1011  MatrixOutline(MatrixOutline::Row),
1012  MatrixCondition())
1013  {
1014  }
1016  explicit RowVector(int n)
1017  : MatrixCommitment
1018  ( MatrixStructure(MatrixStructure::Matrix1d),
1019  MatrixStorage(),
1020  MatrixOutline(MatrixOutline::Row),
1021  MatrixCondition())
1022  {
1023  commitNumCols(n);
1024  }
1025 };
1026 
1028 class MatrixCommitment::Triangular : public MatrixCommitment {
1029 public:
1030  Triangular()
1031  : MatrixCommitment(MatrixStructure::Triangular, MatrixStorage(),
1032  MatrixOutline(), MatrixCondition())
1033  {
1034  }
1035 };
1036 
1039 class MatrixCommitment::Symmetric : public MatrixCommitment {
1040 public:
1041  Symmetric()
1042  : MatrixCommitment(MatrixStructure::Symmetric, MatrixStorage(),
1043  MatrixOutline(), MatrixCondition())
1044  {
1045  }
1046 };
1047 
1051 class MatrixCommitment::Hermitian : public MatrixCommitment {
1052 public:
1053  Hermitian()
1054  : MatrixCommitment
1055  ( MatrixStructure::Hermitian,
1056  MatrixStorage(),
1057  MatrixOutline(),
1058  MatrixCondition().setDiagonal(MatrixCondition::RealDiagonal))
1059  {
1060  }
1061 };
1062 
1066 class MatrixCommitment::SkewSymmetric : public MatrixCommitment {
1067 public:
1068  SkewSymmetric()
1069  : MatrixCommitment
1070  ( MatrixStructure::SkewSymmetric,
1071  MatrixStorage(),
1072  MatrixOutline(),
1073  MatrixCondition().setDiagonal(MatrixCondition::ZeroDiagonal))
1074  {
1075  }
1076 };
1077 
1081 class MatrixCommitment::SkewHermitian : public MatrixCommitment {
1082 public:
1083  SkewHermitian()
1084  : MatrixCommitment
1085  ( MatrixStructure::SkewHermitian,
1086  MatrixStorage(),
1087  MatrixOutline(),
1088  MatrixCondition().setDiagonal(MatrixCondition::ImaginaryDiagonal))
1089  {
1090  }
1091 };
1092 
1095 SimTK_SimTKCOMMON_EXPORT std::ostream&
1096 operator<<(std::ostream& o, const MatrixCommitment&);
1097 
1098 } //namespace SimTK
1099 
1100 #endif // SimTK_SIMMATRIX_MATRIX_CHARACTERISTICS_H_