1 /*
2 *
3 * Template Numerical Toolkit (TNT)
4 *
5 * Mathematical and Computational Sciences Division
6 * National Institute of Technology,
7 * Gaithersburg, MD USA
8 *
9 *
10 * This software was developed at the National Institute of Standards and
11 * Technology (NIST) by employees of the Federal Government in the course
12 * of their official duties. Pursuant to title 17 Section 105 of the
13 * United States Code, this software is not subject to copyright protection
14 * and is in the public domain. NIST assumes no responsibility whatsoever for
15 * its use by other parties, and makes no guarantees, expressed or implied,
16 * about its quality, reliability, or any other characteristic.
17 *
18 */
19
20
21 #ifndef TNT_ARRAY2D_UTILS_H
22 #define TNT_ARRAY2D_UTILS_H
23
24 #include <cstdlib>
25 #include <cassert>
26
27 namespace TNT
28 {
29
30
31 template <class T>
32 std::ostream& operator<<(std::ostream &s, const Array2D<T> &A)
33 {
34 int M=A.dim1();
35 int N=A.dim2();
36
37 s << M << " " << N << "\n";
38
39 for (int i=0; i<M; i++)
40 {
41 for (int j=0; j<N; j++)
42 {
43 s << A[i][j] << " ";
44 }
45 s << "\n";
46 }
47
48
49 return s;
50 }
51
52 template <class T>
53 std::istream& operator>>(std::istream &s, Array2D<T> &A)
54 {
55
56 int M, N;
57
58 s >> M >> N;
59
60 Array2D<T> B(M,N);
61
62 for (int i=0; i<M; i++)
63 for (int j=0; j<N; j++)
64 {
65 s >> B[i][j];
66 }
67
68 A = B;
69 return s;
70 }
71
72
73 template <class T>
74 Array2D<T> operator+(const Array2D<T> &A, const Array2D<T> &B)
75 {
76 int m = A.dim1();
77 int n = A.dim2();
78
79 if (B.dim1() != m || B.dim2() != n )
80 return Array2D<T>();
81
82 else
83 {
84 Array2D<T> C(m,n);
85
86 for (int i=0; i<m; i++)
87 {
88 for (int j=0; j<n; j++)
89 C[i][j] = A[i][j] + B[i][j];
90 }
91 return C;
92 }
93 }
94
95 template <class T>
96 Array2D<T> operator-(const Array2D<T> &A, const Array2D<T> &B)
97 {
98 int m = A.dim1();
99 int n = A.dim2();
100
101 if (B.dim1() != m || B.dim2() != n )
102 return Array2D<T>();
103
104 else
105 {
106 Array2D<T> C(m,n);
107
108 for (int i=0; i<m; i++)
109 {
110 for (int j=0; j<n; j++)
111 C[i][j] = A[i][j] - B[i][j];
112 }
113 return C;
114 }
115 }
116
117
118 template <class T>
119 Array2D<T> operator*(const Array2D<T> &A, const Array2D<T> &B)
120 {
121 int m = A.dim1();
122 int n = A.dim2();
123
124 if (B.dim1() != m || B.dim2() != n )
125 return Array2D<T>();
126
127 else
128 {
129 Array2D<T> C(m,n);
130
131 for (int i=0; i<m; i++)
132 {
133 for (int j=0; j<n; j++)
134 C[i][j] = A[i][j] * B[i][j];
135 }
136 return C;
137 }
138 }
139
140
141
142
143 template <class T>
144 Array2D<T> operator/(const Array2D<T> &A, const Array2D<T> &B)
145 {
146 int m = A.dim1();
147 int n = A.dim2();
148
149 if (B.dim1() != m || B.dim2() != n )
150 return Array2D<T>();
151
152 else
153 {
154 Array2D<T> C(m,n);
155
156 for (int i=0; i<m; i++)
157 {
158 for (int j=0; j<n; j++)
159 C[i][j] = A[i][j] / B[i][j];
160 }
161 return C;
162 }
163 }
164
165
166
167
168
169 template <class T>
170 Array2D<T>& operator+=(Array2D<T> &A, const Array2D<T> &B)
171 {
172 int m = A.dim1();
173 int n = A.dim2();
174
175 if (B.dim1() == m || B.dim2() == n )
176 {
177 for (int i=0; i<m; i++)
178 {
179 for (int j=0; j<n; j++)
180 A[i][j] += B[i][j];
181 }
182 }
183 return A;
184 }
185
186
187
188 template <class T>
189 Array2D<T>& operator-=(Array2D<T> &A, const Array2D<T> &B)
190 {
191 int m = A.dim1();
192 int n = A.dim2();
193
194 if (B.dim1() == m || B.dim2() == n )
195 {
196 for (int i=0; i<m; i++)
197 {
198 for (int j=0; j<n; j++)
199 A[i][j] -= B[i][j];
200 }
201 }
202 return A;
203 }
204
205
206
207 template <class T>
208 Array2D<T>& operator*=(Array2D<T> &A, const Array2D<T> &B)
209 {
210 int m = A.dim1();
211 int n = A.dim2();
212
213 if (B.dim1() == m || B.dim2() == n )
214 {
215 for (int i=0; i<m; i++)
216 {
217 for (int j=0; j<n; j++)
218 A[i][j] *= B[i][j];
219 }
220 }
221 return A;
222 }
223
224
225
226
227
228 template <class T>
229 Array2D<T>& operator/=(Array2D<T> &A, const Array2D<T> &B)
230 {
231 int m = A.dim1();
232 int n = A.dim2();
233
234 if (B.dim1() == m || B.dim2() == n )
235 {
236 for (int i=0; i<m; i++)
237 {
238 for (int j=0; j<n; j++)
239 A[i][j] /= B[i][j];
240 }
241 }
242 return A;
243 }
244
245 /**
246 Matrix Multiply: compute C = A*B, where C[i][j]
247 is the dot-product of row i of A and column j of B.
248
249
250 @param A an (m x n) array
251 @param B an (n x k) array
252 @return the (m x k) array A*B, or a null array (0x0)
253 if the matrices are non-conformant (i.e. the number
254 of columns of A are different than the number of rows of B.)
255
256
257 */
258 template <class T>
matmult(const Array2D<T> & A,const Array2D<T> & B)259 Array2D<T> matmult(const Array2D<T> &A, const Array2D<T> &B)
260 {
261 if (A.dim2() != B.dim1())
262 return Array2D<T>();
263
264 int M = A.dim1();
265 int N = A.dim2();
266 int K = B.dim2();
267
268 Array2D<T> C(M,K);
269
270 for (int i=0; i<M; i++)
271 for (int j=0; j<K; j++)
272 {
273 T sum = 0;
274
275 for (int k=0; k<N; k++)
276 sum += A[i][k] * B [k][j];
277
278 C[i][j] = sum;
279 }
280
281 return C;
282
283 }
284
285 } // namespace TNT
286
287 #endif
288