1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
|
#ifndef GFXMAT2_INCLUDED // -*- C++ -*-
#define GFXMAT2_INCLUDED
#if !defined(__GNUC__)
# pragma once
#endif
/************************************************************************
2x2 Matrix class
$Id: mat2.h 427 2004-09-27 04:45:31Z garland $
************************************************************************/
#include "vec2.h"
namespace gfx
{
class Mat2
{
private:
Vec2 row[2];
public:
// Standard constructors
//
Mat2() { *this = 0.0; }
Mat2(double a, double b, double c, double d)
{ row[0][0]=a; row[0][1]=b; row[1][0]=c; row[1][1]=d; }
Mat2(const Vec2 &r0,const Vec2 &r1) { row[0]=r0; row[1]=r1; }
Mat2(const Mat2 &m) { *this = m; }
// Descriptive interface
//
typedef double value_type;
typedef Vec2 vector_type;
typedef Mat2 inverse_type;
static int dim() { return 2; }
// Access methods note: A(i, j) == row i, col j
//
double& operator()(int i, int j) { return row[i][j]; }
double operator()(int i, int j) const { return row[i][j]; }
Vec2& operator[](int i) { return row[i]; }
const Vec2& operator[](int i) const { return row[i]; }
inline Vec2 col(int i) const {return Vec2(row[0][i],row[1][i]);}
operator double*() { return row[0]; }
operator const double*() { return row[0]; }
operator const double*() const { return row[0]; }
// Assignment methods
//
inline Mat2& operator=(const Mat2& m);
inline Mat2& operator=(double s);
inline Mat2& operator+=(const Mat2& m);
inline Mat2& operator-=(const Mat2& m);
inline Mat2& operator*=(double s);
inline Mat2& operator/=(double s);
// Construction of standard matrices
//
static Mat2 I();
static Mat2 outer_product(const Vec2 &u, const Vec2 &v)
{ return Mat2(u[0]*v[0], u[0]*v[1], u[1]*v[0], u[1]*v[1]); }
static Mat2 outer_product(const Vec2 &u) { return outer_product(u,u); }
Mat2 &diag(double d);
Mat2 &ident() { return diag(1.0); }
};
////////////////////////////////////////////////////////////////////////
//
// Method definitions
//
inline Mat2& Mat2::operator=(const Mat2& m)
{ row[0]=m[0]; row[1]=m[1]; return *this; }
inline Mat2& Mat2::operator=(double s)
{ row[0]=s; row[1]=s; return *this; }
inline Mat2& Mat2::operator+=(const Mat2& m)
{ row[0] += m.row[0]; row[1] += m.row[1]; return *this;}
inline Mat2& Mat2::operator-=(const Mat2& m)
{ row[0] -= m.row[0]; row[1] -= m.row[1]; return *this; }
inline Mat2& Mat2::operator*=(double s)
{ row[0] *= s; row[1] *= s; return *this; }
inline Mat2& Mat2::operator/=(double s)
{ row[0] /= s; row[1] /= s; return *this; }
////////////////////////////////////////////////////////////////////////
//
// Operator definitions
//
inline Mat2 operator+(const Mat2 &n, const Mat2 &m)
{ return Mat2(n[0]+m[0], n[1]+m[1]); }
inline Mat2 operator-(const Mat2 &n, const Mat2 &m)
{ return Mat2(n[0]-m[0], n[1]-m[1]); }
inline Mat2 operator-(const Mat2 &m)
{ return Mat2(-m[0], -m[1]); }
inline Mat2 operator*(double s, const Mat2 &m)
{ return Mat2(m[0]*s, m[1]*s); }
inline Mat2 operator*(const Mat2 &m, double s)
{ return s*m; }
inline Mat2 operator/(const Mat2 &m, double s)
{ return Mat2(m[0]/s, m[1]/s); }
inline Vec2 operator*(const Mat2 &m, const Vec2 &v)
{ return Vec2(m[0]*v, m[1]*v); }
extern Mat2 operator*(const Mat2 &n, const Mat2 &m);
inline std::ostream &operator<<(std::ostream &out, const Mat2& M)
{ return out << M[0] << std::endl << M[1]; }
inline std::istream &operator>>(std::istream &in, Mat2& M)
{ return in >> M[0] >> M[1]; }
////////////////////////////////////////////////////////////////////////
//
// Misc. function definitions
//
inline double det(const Mat2 &m)
{ return m(0,0)*m(1,1) - m(0,1)*m(1,0); }
inline double trace(const Mat2 &m)
{ return m(0,0) + m(1,1); }
inline Mat2 transpose(const Mat2 &m)
{ return Mat2(m.col(0), m.col(1)); }
inline Mat2 adjoint(const Mat2 &m)
{ return Mat2(perp(m[1]), -perp(m[0])); }
extern double invert(Mat2 &m_inv, const Mat2 &m);
extern bool eigenvalues(const Mat2&, Vec2& evals);
extern bool eigenvectors(const Mat2&, const Vec2& evals, Vec2 evecs[2]);
extern bool eigen(const Mat2&, Vec2& evals, Vec2 evecs[2]);
} // namespace gfx
// GFXMAT2_INCLUDED
#endif
|
