111
|
1 // { dg-do run }
|
|
2 // { dg-options "-fpic" { target fpic } }
|
|
3
|
|
4 typedef __SIZE_TYPE__ size_t;
|
|
5
|
|
6 template <typename NumType>
|
|
7 inline
|
|
8 NumType
|
|
9 absolute(NumType const& x)
|
|
10 {
|
|
11 if (x < NumType(0)) return -x;
|
|
12 return x;
|
|
13 }
|
|
14
|
|
15 class trivial_accessor
|
|
16 {
|
|
17 public:
|
|
18 typedef size_t index_type;
|
|
19 struct index_value_type {};
|
|
20
|
|
21 trivial_accessor() : size_(0) {}
|
|
22
|
|
23 trivial_accessor(size_t const& n) : size_(n) {}
|
|
24
|
|
25 size_t size_1d() const { return size_; }
|
|
26
|
|
27 protected:
|
|
28 size_t size_;
|
|
29 };
|
|
30
|
|
31 namespace N0
|
|
32 {
|
|
33 template <typename ElementType,
|
|
34 typename AccessorType = trivial_accessor>
|
|
35 class const_ref
|
|
36 {
|
|
37 public:
|
|
38 typedef ElementType value_type;
|
|
39 typedef size_t size_type;
|
|
40
|
|
41 typedef AccessorType accessor_type;
|
|
42 typedef typename accessor_type::index_type index_type;
|
|
43 typedef typename accessor_type::index_value_type index_value_type;
|
|
44
|
|
45 const_ref() {}
|
|
46
|
|
47 const_ref(const ElementType* begin, accessor_type const& accessor)
|
|
48 : begin_(begin), accessor_(accessor)
|
|
49 {
|
|
50 init();
|
|
51 }
|
|
52
|
|
53 const_ref(const ElementType* begin, index_value_type const& n0)
|
|
54 : begin_(begin), accessor_(n0)
|
|
55 {
|
|
56 init();
|
|
57 }
|
|
58
|
|
59 const_ref(const ElementType* begin, index_value_type const& n0,
|
|
60 index_value_type const& n1)
|
|
61 : begin_(begin), accessor_(n0, n1)
|
|
62 {
|
|
63 init();
|
|
64 }
|
|
65
|
|
66 const_ref(const ElementType* begin, index_value_type const& n0,
|
|
67 index_value_type const& n1,
|
|
68 index_value_type const& n2)
|
|
69 : begin_(begin), accessor_(n0, n1, n2)
|
|
70 {
|
|
71 init();
|
|
72 }
|
|
73
|
|
74 accessor_type const& accessor() const { return accessor_; }
|
|
75 size_type size() const { return size_; }
|
|
76
|
|
77 const ElementType* begin() const { return begin_; }
|
|
78 const ElementType* end() const { return end_; }
|
|
79
|
|
80 ElementType const&
|
|
81 operator[](size_type i) const { return begin_[i]; }
|
|
82
|
|
83 const_ref<ElementType>
|
|
84 as_1d() const
|
|
85 {
|
|
86 return const_ref<ElementType>(begin_, size_);
|
|
87 }
|
|
88
|
|
89 protected:
|
|
90 void
|
|
91 init()
|
|
92 {
|
|
93 size_ = accessor_.size_1d();
|
|
94 end_ = begin_ + size_;
|
|
95 }
|
|
96
|
|
97 const ElementType* begin_;
|
|
98 accessor_type accessor_;
|
|
99 size_type size_;
|
|
100 const ElementType* end_;
|
|
101 };
|
|
102 }
|
|
103
|
|
104 template <typename ElementType,
|
|
105 typename AccessorType = trivial_accessor>
|
|
106 class ref : public N0::const_ref<ElementType, AccessorType>
|
|
107 {
|
|
108 public:
|
|
109 typedef ElementType value_type;
|
|
110 typedef size_t size_type;
|
|
111
|
|
112 typedef N0::const_ref<ElementType, AccessorType> base_class;
|
|
113 typedef AccessorType accessor_type;
|
|
114 typedef typename accessor_type::index_type index_type;
|
|
115
|
|
116 ref() {}
|
|
117
|
|
118 ElementType*
|
|
119 begin() const { return const_cast<ElementType*>(this->begin_); }
|
|
120
|
|
121 ElementType*
|
|
122 end() const { return const_cast<ElementType*>(this->end_); }
|
|
123
|
|
124 ElementType&
|
|
125 operator[](size_type i) const { return begin()[i]; }
|
|
126 };
|
|
127
|
|
128 namespace N1 {
|
|
129 template <typename ElementType, size_t N>
|
|
130 class tiny_plain
|
|
131 {
|
|
132 public:
|
|
133 typedef ElementType value_type;
|
|
134 typedef size_t size_type;
|
|
135
|
|
136 static const size_t fixed_size=N;
|
|
137
|
|
138 ElementType elems[N];
|
|
139
|
|
140 tiny_plain() {}
|
|
141
|
|
142 static size_type size() { return N; }
|
|
143
|
|
144 ElementType* begin() { return elems; }
|
|
145 const ElementType* begin() const { return elems; }
|
|
146 ElementType* end() { return elems+N; }
|
|
147 const ElementType* end() const { return elems+N; }
|
|
148 ElementType& operator[](size_type i) { return elems[i]; }
|
|
149 ElementType const& operator[](size_type i) const { return elems[i]; }
|
|
150 };
|
|
151
|
|
152 template <typename ElementType, size_t N>
|
|
153 class tiny : public tiny_plain<ElementType, N>
|
|
154 {
|
|
155 public:
|
|
156 typedef ElementType value_type;
|
|
157 typedef size_t size_type;
|
|
158
|
|
159 typedef tiny_plain<ElementType, N> base_class;
|
|
160
|
|
161 tiny() {}
|
|
162 };
|
|
163 }
|
|
164
|
|
165 template <typename NumType>
|
|
166 class mat3 : public N1::tiny_plain<NumType, 9>
|
|
167 {
|
|
168 public:
|
|
169 typedef typename N1::tiny_plain<NumType, 9> base_type;
|
|
170
|
|
171 mat3() {}
|
|
172 mat3(NumType const& e00, NumType const& e01, NumType const& e02,
|
|
173 NumType const& e10, NumType const& e11, NumType const& e12,
|
|
174 NumType const& e20, NumType const& e21, NumType const& e22)
|
|
175 : base_type(e00, e01, e02, e10, e11, e12, e20, e21, e22)
|
|
176 {}
|
|
177 mat3(base_type const& a)
|
|
178 : base_type(a)
|
|
179 {}
|
|
180
|
|
181 NumType const&
|
|
182 operator()(size_t r, size_t c) const
|
|
183 {
|
|
184 return this->elems[r * 3 + c];
|
|
185 }
|
|
186 NumType&
|
|
187 operator()(size_t r, size_t c)
|
|
188 {
|
|
189 return this->elems[r * 3 + c];
|
|
190 }
|
|
191
|
|
192 NumType
|
|
193 trace() const
|
|
194 {
|
|
195 mat3 const& m = *this;
|
|
196 return m[0] + m[4] + m[8];
|
|
197 }
|
|
198
|
|
199 NumType
|
|
200 determinant() const
|
|
201 {
|
|
202 mat3 const& m = *this;
|
|
203 return m[0] * (m[4] * m[8] - m[5] * m[7])
|
|
204 - m[1] * (m[3] * m[8] - m[5] * m[6])
|
|
205 + m[2] * (m[3] * m[7] - m[4] * m[6]);
|
|
206 }
|
|
207 };
|
|
208
|
|
209 template <typename NumType>
|
|
210 inline
|
|
211 mat3<NumType>
|
|
212 operator-(mat3<NumType> const& v)
|
|
213 {
|
|
214 mat3<NumType> result;
|
|
215 for(size_t i=0;i<9;i++) {
|
|
216 result[i] = -v[i];
|
|
217 }
|
|
218 return result;
|
|
219 }
|
|
220
|
|
221 class mat_grid : public N1::tiny<size_t, 2>
|
|
222 {
|
|
223 public:
|
|
224 typedef N1::tiny<size_t, 2> index_type;
|
|
225 typedef index_type::value_type index_value_type;
|
|
226
|
|
227 mat_grid() { this->elems[0]=0; this->elems[1]=0; }
|
|
228
|
|
229 mat_grid(index_type const& n) : index_type(n) {}
|
|
230
|
|
231 mat_grid(index_value_type const& n0, index_value_type const& n1)
|
|
232 { this->elems[0]=n0; this->elems[1]=n1; }
|
|
233
|
|
234 size_t size_1d() const { return elems[0] * elems[1]; }
|
|
235
|
|
236 size_t
|
|
237 operator()(index_value_type const& r, index_value_type const& c) const
|
|
238 {
|
|
239 return r * elems[1] + c;
|
|
240 }
|
|
241 };
|
|
242
|
|
243 template <typename NumType, typename AccessorType = mat_grid>
|
|
244 class mat_const_ref : public N0::const_ref<NumType, AccessorType>
|
|
245 {
|
|
246 public:
|
|
247 typedef AccessorType accessor_type;
|
|
248 typedef typename N0::const_ref<NumType, AccessorType> base_type;
|
|
249 typedef typename accessor_type::index_value_type index_value_type;
|
|
250
|
|
251 mat_const_ref() {}
|
|
252
|
|
253 mat_const_ref(const NumType* begin, accessor_type const& grid)
|
|
254 : base_type(begin, grid)
|
|
255 {}
|
|
256
|
|
257 mat_const_ref(const NumType* begin, index_value_type const& n_rows,
|
|
258 index_value_type const& n_columns)
|
|
259 : base_type(begin, accessor_type(n_rows, n_columns))
|
|
260 {}
|
|
261
|
|
262 accessor_type
|
|
263 grid() const { return this->accessor(); }
|
|
264
|
|
265 index_value_type const&
|
|
266 n_rows() const { return this->accessor()[0]; }
|
|
267
|
|
268 index_value_type const&
|
|
269 n_columns() const { return this->accessor()[1]; }
|
|
270
|
|
271 NumType const&
|
|
272 operator()(index_value_type const& r, index_value_type const& c) const
|
|
273 {
|
|
274 return this->begin()[this->accessor()(r, c)];
|
|
275 }
|
|
276 };
|
|
277
|
|
278 template <typename NumType, typename AccessorType = mat_grid>
|
|
279 class mat_ref : public mat_const_ref<NumType, AccessorType>
|
|
280 {
|
|
281 public:
|
|
282 typedef AccessorType accessor_type;
|
|
283 typedef mat_const_ref<NumType, AccessorType> base_type;
|
|
284 typedef typename accessor_type::index_value_type index_value_type;
|
|
285
|
|
286 mat_ref() {}
|
|
287
|
|
288 mat_ref(NumType* begin, accessor_type const& grid)
|
|
289 : base_type(begin, grid)
|
|
290 {}
|
|
291
|
|
292 mat_ref(NumType* begin, index_value_type n_rows,
|
|
293 index_value_type n_columns)
|
|
294 : base_type(begin, accessor_type(n_rows, n_columns))
|
|
295 {}
|
|
296
|
|
297 NumType*
|
|
298 begin() const { return const_cast<NumType*>(this->begin_); }
|
|
299
|
|
300 NumType*
|
|
301 end() const { return const_cast<NumType*>(this->end_); }
|
|
302
|
|
303 NumType&
|
|
304 operator[](index_value_type const& i) const { return begin()[i]; }
|
|
305
|
|
306 NumType&
|
|
307 operator()(index_value_type const& r, index_value_type const& c) const
|
|
308 {
|
|
309 return this->begin()[this->accessor()(r, c)];
|
|
310 }
|
|
311 };
|
|
312
|
|
313 template <typename AnyType>
|
|
314 inline void
|
|
315 swap(AnyType* a, AnyType* b, size_t n)
|
|
316 {
|
|
317 for(size_t i=0;i<n;i++) {
|
|
318 AnyType t = a[i]; a[i] = b[i]; b[i] = t;
|
|
319 }
|
|
320 }
|
|
321
|
|
322 template <typename IntType>
|
|
323 size_t
|
|
324 form_t(mat_ref<IntType>& m,
|
|
325 mat_ref<IntType> const& t)
|
|
326 {
|
|
327 typedef size_t size_t;
|
|
328 size_t mr = m.n_rows();
|
|
329 size_t mc = m.n_columns();
|
|
330 size_t tc = t.n_columns();
|
|
331 if (tc) {
|
|
332 }
|
|
333 size_t i, j;
|
|
334 for (i = j = 0; i < mr && j < mc;) {
|
|
335 size_t k = i; while (k < mr && m(k,j) == 0) k++;
|
|
336 if (k == mr)
|
|
337 j++;
|
|
338 else {
|
|
339 if (i != k) {
|
|
340 swap(&m(i,0), &m(k,0), mc);
|
|
341 if (tc) swap(&t(i,0), &t(k,0), tc);
|
|
342 }
|
|
343 for (k++; k < mr; k++) {
|
|
344 IntType a = absolute(m(k, j));
|
|
345 if (a != 0 && a < absolute(m(i,j))) {
|
|
346 swap(&m(i,0), &m(k,0), mc);
|
|
347 if (tc) swap(&t(i,0), &t(k,0), tc);
|
|
348 }
|
|
349 }
|
|
350 if (m(i,j) < 0) {
|
|
351 for(size_t ic=0;ic<mc;ic++) m(i,ic) *= -1;
|
|
352 if (tc) for(size_t ic=0;ic<tc;ic++) t(i,ic) *= -1;
|
|
353 }
|
|
354 bool cleared = true;
|
|
355 for (k = i+1; k < mr; k++) {
|
|
356 IntType a = m(k,j) / m(i,j);
|
|
357 if (a != 0) {
|
|
358 for(size_t ic=0;ic<mc;ic++) m(k,ic) -= a * m(i,ic);
|
|
359 if (tc) for(size_t ic=0;ic<tc;ic++) t(k,ic) -= a * t(i,ic);
|
|
360 }
|
|
361 if (m(k,j) != 0) cleared = false;
|
|
362 }
|
|
363 if (cleared) { i++; j++; }
|
|
364 }
|
|
365 }
|
|
366 m = mat_ref<IntType>(m.begin(), i, mc);
|
|
367 return i;
|
|
368 }
|
|
369
|
|
370 template <typename IntType>
|
|
371 size_t
|
|
372 form(mat_ref<IntType>& m)
|
|
373 {
|
|
374 mat_ref<IntType> t(0,0,0);
|
|
375 return form_t(m, t);
|
|
376 }
|
|
377
|
|
378 typedef mat3<int> sg_mat3;
|
|
379
|
|
380 class rot_mx
|
|
381 {
|
|
382 public:
|
|
383 explicit
|
|
384 rot_mx(sg_mat3 const& m, int denominator=1)
|
|
385 : num_(m), den_(denominator)
|
|
386 {}
|
|
387
|
|
388 sg_mat3 const&
|
|
389 num() const { return num_; }
|
|
390 sg_mat3&
|
|
391 num() { return num_; }
|
|
392
|
|
393 int const&
|
|
394 operator[](size_t i) const { return num_[i]; }
|
|
395 int&
|
|
396 operator[](size_t i) { return num_[i]; }
|
|
397
|
|
398 int
|
|
399 const& operator()(int r, int c) const { return num_(r, c); }
|
|
400 int&
|
|
401 operator()(int r, int c) { return num_(r, c); }
|
|
402
|
|
403 int const&
|
|
404 den() const { return den_; }
|
|
405 int&
|
|
406 den() { return den_; }
|
|
407
|
|
408 rot_mx
|
|
409 minus_unit_mx() const
|
|
410 {
|
|
411 rot_mx result(*this);
|
|
412 for (size_t i=0;i<9;i+=4) result[i] -= den_;
|
|
413 return result;
|
|
414 }
|
|
415
|
|
416 rot_mx
|
|
417 operator-() const { return rot_mx(-num_, den_); }
|
|
418
|
|
419 int
|
|
420 type() const;
|
|
421
|
|
422 int
|
|
423 order(int type=0) const;
|
|
424
|
|
425 private:
|
|
426 sg_mat3 num_;
|
|
427 int den_;
|
|
428 };
|
|
429
|
|
430 class rot_mx_info
|
|
431 {
|
|
432 public:
|
|
433 rot_mx_info(rot_mx const& r);
|
|
434
|
|
435 int type() const { return type_; }
|
|
436
|
|
437 private:
|
|
438 int type_;
|
|
439 };
|
|
440
|
|
441 int rot_mx::type() const
|
|
442 {
|
|
443 int det = num_.determinant();
|
|
444 if (det == -1 || det == 1) {
|
|
445 switch (num_.trace()) {
|
|
446 case -3: return -1;
|
|
447 case -2: return -6;
|
|
448 case -1: if (det == -1) return -4;
|
|
449 else return 2;
|
|
450 case 0: if (det == -1) return -3;
|
|
451 else return 3;
|
|
452 case 1: if (det == -1) return -2;
|
|
453 else return 4;
|
|
454 case 2: return 6;
|
|
455 case 3: return 1;
|
|
456 }
|
|
457 }
|
|
458 return 0;
|
|
459 }
|
|
460
|
|
461 int rot_mx::order(int type) const
|
|
462 {
|
|
463 if (type == 0) type = rot_mx::type();
|
|
464 if (type > 0) return type;
|
|
465 if (type % 2) return -type * 2;
|
|
466 return -type;
|
|
467 }
|
|
468
|
|
469 rot_mx_info::rot_mx_info(rot_mx const& r)
|
|
470 : type_(r.type())
|
|
471 {
|
|
472 if (type_ == 0) {
|
|
473 return;
|
|
474 }
|
|
475 rot_mx proper_r = r;
|
|
476 int proper_order = type_;
|
|
477 // THE PROBLEM IS AROUND HERE
|
|
478 if (proper_order < 0) {
|
|
479 proper_order *= -1;
|
|
480 proper_r = -proper_r; // THIS FAILS ...
|
|
481 }
|
|
482 if (proper_order > 1) {
|
|
483 rot_mx rmi = proper_r.minus_unit_mx(); // ... THEREFORE WRONG HERE
|
|
484 mat_ref<int> re_mx(rmi.num().begin(), 3, 3);
|
|
485 if (form(re_mx) != 2) {
|
|
486 type_ = 0;
|
|
487 }
|
|
488 }
|
|
489 }
|
|
490
|
|
491 int main()
|
|
492 {
|
|
493 N1::tiny<int, 9> e;
|
|
494 e[0] = 1; e[1] = 0; e[2] = 0;
|
|
495 e[3] = 0; e[4] = -1; e[5] = 0;
|
|
496 e[6] = 0; e[7] = 0; e[8] = 1;
|
|
497 rot_mx r(e);
|
|
498 rot_mx_info ri(r);
|
|
499 if (ri.type() != -2)
|
|
500 __builtin_abort ();
|
|
501 return 0;
|
|
502 }
|