Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

hermitian.hpp @ 857

Revision 857, 93.7 KB checked in by igarcia, 19 years ago (diff)

Line
1	//
2	// Copyright (c) 2000-2002
3	// Joerg Walter, Mathias Koch
4	//
5	// Permission to use, copy, modify, distribute and sell this software
6	// and its documentation for any purpose is hereby granted without fee,
7	// provided that the above copyright notice appear in all copies and
8	// that both that copyright notice and this permission notice appear
9	// in supporting documentation. The authors make no representations
10	// about the suitability of this software for any purpose.
11	// It is provided "as is" without express or implied warranty.
12	//
13	// The authors gratefully acknowledge the support of
14	// GeNeSys mbH & Co. KG in producing this work.
15	//
16
17	#ifndef BOOST_UBLAS_HERMITIAN_H
18	#define BOOST_UBLAS_HERMITIAN_H
19
20	#include <boost/numeric/ublas/matrix.hpp>
21	#include <boost/numeric/ublas/detail/temporary.hpp>
22
23	// Iterators based on ideas of Jeremy Siek
24	// Hermitian matrices are square. Thanks to Peter Schmitteckert for spotting this.
25
26	namespace boost { namespace numeric { namespace ublas {
27
28	template<class M>
29	bool is_hermitian (const M &m) {
30	typedef typename M::size_type size_type;
31
32	if (m.size1 () != m.size2 ())
33	return false;
34	size_type size = BOOST_UBLAS_SAME (m.size1 (), m.size2 ());
35	for (size_type i = 0; i < size; ++ i) {
36	for (size_type j = i; j < size; ++ j) {
37	if (m (i, j) != conj (m (j, i)))
38	return false;
39	}
40	}
41	return true;
42	}
43
44	#ifdef BOOST_UBLAS_STRICT_HERMITIAN
45
46	template<class M>
47	class hermitian_matrix_element:
48	public container_reference<M> {
49	public:
50	typedef M matrix_type;
51	typedef typename M::size_type size_type;
52	typedef typename M::value_type value_type;
53	typedef const value_type &const_reference;
54	typedef value_type &reference;
55	typedef value_type *pointer;
56
57	// Construction and destruction
58	BOOST_UBLAS_INLINE
59	hermitian_matrix_element (matrix_type &m, size_type i, size_type j, value_type d):
60	container_reference<matrix_type> (m), i_ (i), j_ (j), d_ (d), dirty_ (false) {}
61	BOOST_UBLAS_INLINE
62	hermitian_matrix_element (const hermitian_matrix_element &p):
63	container_reference<matrix_type> (p), i_ (p.i_), d_ (p.d_), dirty_ (p.dirty_) {}
64	BOOST_UBLAS_INLINE
65	~hermitian_matrix_element () {
66	if (dirty_)
67	(*this) ().insert_element (i_, j_, d_);
68	}
69
70	// Assignment
71	BOOST_UBLAS_INLINE
72	hermitian_matrix_element &operator = (const hermitian_matrix_element &p) {
73	// Overide the implict copy assignment
74	d_ = p.d_;
75	dirty_ = true;
76	return *this;
77	}
78	template<class D>
79	BOOST_UBLAS_INLINE
80	hermitian_matrix_element &operator = (const D &d) {
81	d_ = d;
82	dirty_ = true;
83	return *this;
84	}
85	template<class D>
86	BOOST_UBLAS_INLINE
87	hermitian_matrix_element &operator += (const D &d) {
88	d_ += d;
89	dirty_ = true;
90	return *this;
91	}
92	template<class D>
93	BOOST_UBLAS_INLINE
94	hermitian_matrix_element &operator -= (const D &d) {
95	d_ -= d;
96	dirty_ = true;
97	return *this;
98	}
99	template<class D>
100	BOOST_UBLAS_INLINE
101	hermitian_matrix_element &operator *= (const D &d) {
102	d_ *= d;
103	dirty_ = true;
104	return *this;
105	}
106	template<class D>
107	BOOST_UBLAS_INLINE
108	hermitian_matrix_element &operator /= (const D &d) {
109	d_ /= d;
110	dirty_ = true;
111	return *this;
112	}
113
114	// Comparison
115	template<class D>
116	BOOST_UBLAS_INLINE
117	bool operator == (const D &d) const {
118	return d_ == d;
119	}
120	template<class D>
121	BOOST_UBLAS_INLINE
122	bool operator != (const D &d) const {
123	return d_ != d;
124	}
125
126	// Conversion
127	BOOST_UBLAS_INLINE
128	operator const_reference () const {
129	return d_;
130	}
131
132	// Swapping
133	BOOST_UBLAS_INLINE
134	void swap (hermitian_matrix_element p) {
135	if (this != &p) {
136	dirty_ = true;
137	p.dirty_ = true;
138	std::swap (d_, p.d_);
139	}
140	}
141	BOOST_UBLAS_INLINE
142	friend void swap (hermitian_matrix_element p1, hermitian_matrix_element p2) {
143	p1.swap (p2);
144	}
145
146	private:
147	size_type i_;
148	size_type j_;
149	value_type d_;
150	bool dirty_;
151	};
152
153	template<class M>
154	struct type_traits<hermitian_matrix_element<M> > {
155	typedef typename M::value_type element_type;
156	typedef type_traits<hermitian_matrix_element<M> > self_type;
157	typedef typename type_traits<element_type>::value_type value_type;
158	typedef typename type_traits<element_type>::const_reference const_reference;
159	typedef hermitian_matrix_element<M> reference;
160	typedef typename type_traits<element_type>::real_type real_type;
161	typedef typename type_traits<element_type>::precision_type precision_type;
162
163	static const unsigned plus_complexity = type_traits<element_type>::plus_complexity;
164	static const unsigned multiplies_complexity = type_traits<element_type>::multiplies_complexity;
165
166	static
167	BOOST_UBLAS_INLINE
168	real_type real (const_reference t) {
169	return type_traits<element_type>::real (t);
170	}
171	static
172	BOOST_UBLAS_INLINE
173	real_type imag (const_reference t) {
174	return type_traits<element_type>::imag (t);
175	}
176	static
177	BOOST_UBLAS_INLINE
178	value_type conj (const_reference t) {
179	return type_traits<element_type>::conj (t);
180	}
181
182	static
183	BOOST_UBLAS_INLINE
184	real_type abs (const_reference t) {
185	return type_traits<element_type>::abs (t);
186	}
187	static
188	BOOST_UBLAS_INLINE
189	value_type sqrt (const_reference t) {
190	return type_traits<element_type>::sqrt (t);
191	}
192
193	static
194	BOOST_UBLAS_INLINE
195	real_type norm_1 (const_reference t) {
196	return type_traits<element_type>::norm_1 (t);
197	}
198	static
199	BOOST_UBLAS_INLINE
200	real_type norm_2 (const_reference t) {
201	return type_traits<element_type>::norm_2 (t);
202	}
203	static
204	BOOST_UBLAS_INLINE
205	real_type norm_inf (const_reference t) {
206	return type_traits<element_type>::norm_inf (t);
207	}
208
209	static
210	BOOST_UBLAS_INLINE
211	bool equals (const_reference t1, const_reference t2) {
212	return type_traits<element_type>::equals (t1, t2);
213	}
214	};
215
216	template<class M1, class T2>
217	struct promote_traits<hermitian_matrix_element<M1>, T2> {
218	typedef typename promote_traits<typename hermitian_matrix_element<M1>::value_type, T2>::promote_type promote_type;
219	};
220	template<class T1, class M2>
221	struct promote_traits<T1, hermitian_matrix_element<M2> > {
222	typedef typename promote_traits<T1, typename hermitian_matrix_element<M2>::value_type>::promote_type promote_type;
223	};
224	template<class M1, class M2>
225	struct promote_traits<hermitian_matrix_element<M1>, hermitian_matrix_element<M2> > {
226	typedef typename promote_traits<typename hermitian_matrix_element<M1>::value_type,
227	typename hermitian_matrix_element<M2>::value_type>::promote_type promote_type;
228	};
229
230	#endif
231
232	// Array based hermitian matrix class
233	template<class T, class TRI, class L, class A>
234	class hermitian_matrix:
235	public matrix_container<hermitian_matrix<T, TRI, L, A> > {
236
237	typedef T &true_reference;
238	typedef T *pointer;
239	typedef TRI triangular_type;
240	typedef L layout_type;
241	typedef hermitian_matrix<T, TRI, L, A> self_type;
242	public:
243	#ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
244	using matrix_container<self_type>::operator ();
245	#endif
246	typedef typename A::size_type size_type;
247	typedef typename A::difference_type difference_type;
248	typedef T value_type;
249	// FIXME no better way to not return the address of a temporary?
250	// typedef const T &const_reference;
251	typedef const T const_reference;
252	#ifndef BOOST_UBLAS_STRICT_HERMITIAN
253	typedef T &reference;
254	#else
255	typedef hermitian_matrix_element<self_type> reference;
256	#endif
257	typedef A array_type;
258
259	typedef const matrix_reference<const self_type> const_closure_type;
260	typedef matrix_reference<self_type> closure_type;
261	typedef vector<T, A> vector_temporary_type;
262	typedef matrix<T, L, A> matrix_temporary_type; // general sub-matrix
263	typedef packed_tag storage_category;
264	typedef typename L::orientation_category orientation_category;
265
266	// Construction and destruction
267	BOOST_UBLAS_INLINE
268	hermitian_matrix ():
269	matrix_container<self_type> (),
270	size_ (0), data_ (0) {}
271	BOOST_UBLAS_INLINE
272	hermitian_matrix (size_type size):
273	matrix_container<self_type> (),
274	size_ (BOOST_UBLAS_SAME (size, size)), data_ (triangular_type::packed_size (layout_type (), size, size)) {
275	}
276	BOOST_UBLAS_INLINE
277	hermitian_matrix (size_type size1, size_type size2):
278	matrix_container<self_type> (),
279	size_ (BOOST_UBLAS_SAME (size1, size2)), data_ (triangular_type::packed_size (layout_type (), size1, size2)) {
280	}
281	BOOST_UBLAS_INLINE
282	hermitian_matrix (size_type size, const array_type &data):
283	matrix_container<self_type> (),
284	size_ (size), data_ (data) {}
285	BOOST_UBLAS_INLINE
286	hermitian_matrix (const hermitian_matrix &m):
287	matrix_container<self_type> (),
288	size_ (m.size_), data_ (m.data_) {}
289	template<class AE>
290	BOOST_UBLAS_INLINE
291	hermitian_matrix (const matrix_expression<AE> &ae):
292	matrix_container<self_type> (),
293	size_ (BOOST_UBLAS_SAME (ae ().size1 (), ae ().size2 ())),
294	data_ (triangular_type::packed_size (layout_type (), size_, size_)) {
295	matrix_assign<scalar_assign> (*this, ae);
296	}
297
298	// Accessors
299	BOOST_UBLAS_INLINE
300	size_type size1 () const {
301	return size_;
302	}
303	BOOST_UBLAS_INLINE
304	size_type size2 () const {
305	return size_;
306	}
307
308	// Storage accessors
309	BOOST_UBLAS_INLINE
310	const array_type &data () const {
311	return data_;
312	}
313	BOOST_UBLAS_INLINE
314	array_type &data () {
315	return data_;
316	}
317
318	// Resizing
319	BOOST_UBLAS_INLINE
320	void resize (size_type size, bool preserve = true) {
321	if (preserve) {
322	self_type temporary (size, size);
323	detail::matrix_resize_preserve<layout_type> (*this, temporary);
324	}
325	else {
326	data ().resize (triangular_type::packed_size (layout_type (), size, size));
327	size_ = size;
328	}
329	}
330	BOOST_UBLAS_INLINE
331	void resize (size_type size1, size_type size2, bool preserve = true) {
332	resize (BOOST_UBLAS_SAME (size1, size2), preserve);
333	}
334	BOOST_UBLAS_INLINE
335	void resize_packed_preserve (size_type size) {
336	size_ = BOOST_UBLAS_SAME (size, size);
337	data ().resize (triangular_type::packed_size (layout_type (), size_, size_), value_type ());
338	}
339
340	// Element access
341	BOOST_UBLAS_INLINE
342	const_reference operator () (size_type i, size_type j) const {
343	BOOST_UBLAS_CHECK (i < size_, bad_index ());
344	BOOST_UBLAS_CHECK (j < size_, bad_index ());
345	// if (i == j)
346	// return type_traits<value_type>::real (data () [triangular_type::element (layout_type (), i, size_, i, size_)]);
347	// else
348	if (triangular_type::other (i, j))
349	return data () [triangular_type::element (layout_type (), i, size_, j, size_)];
350	else
351	return type_traits<value_type>::conj (data () [triangular_type::element (layout_type (), j, size_, i, size_)]);
352	}
353	BOOST_UBLAS_INLINE
354	true_reference at_element (size_type i, size_type j) {
355	BOOST_UBLAS_CHECK (i < size_, bad_index ());
356	BOOST_UBLAS_CHECK (j < size_, bad_index ());
357	BOOST_UBLAS_CHECK (triangular_type::other (i, j), bad_index ());
358	return data () [triangular_type::element (layout_type (), i, size_, j, size_)];
359	}
360	BOOST_UBLAS_INLINE
361	reference operator () (size_type i, size_type j) {
362	#ifndef BOOST_UBLAS_STRICT_HERMITIAN
363	if (triangular_type::other (i, j))
364	return at_element (i, j);
365	else {
366	external_logic ().raise ();
367	// arbitary return value
368	return data () [triangular_type::element (layout_type (), j, size_, i, size_)];
369	}
370	#else
371	if (triangular_type::other (i, j))
372	return reference (*this, i, j, data () [triangular_type::element (layout_type (), i, size_, j, size_)]);
373	else
374	return reference (*this, i, j, type_traits<value_type>::conj (data () [triangular_type::element (layout_type (), j, size_, i, size_)]));
375	#endif
376	}
377
378	// Element assignemnt
379	BOOST_UBLAS_INLINE
380	true_reference insert_element (size_type i, size_type j, const_reference t) {
381	BOOST_UBLAS_CHECK (i < size_, bad_index ());
382	BOOST_UBLAS_CHECK (j < size_, bad_index ());
383	if (triangular_type::other (i, j)) {
384	return (data () [triangular_type::element (layout_type (), i, size_, j, size_)] = t);
385	} else {
386	return (data () [triangular_type::element (layout_type (), j, size_, i, size_)] = type_traits<value_type>::conj (t));
387	}
388	}
389	BOOST_UBLAS_INLINE
390	void erase_element (size_type i, size_type j) {
391	BOOST_UBLAS_CHECK (i < size_, bad_index ());
392	BOOST_UBLAS_CHECK (j < size_, bad_index ());
393	data () [triangular_type::element (layout_type (), i, size_, j, size_)] = value_type/zero/();
394	}
395
396	// Zeroing
397	BOOST_UBLAS_INLINE
398	void clear () {
399	std::fill (data ().begin (), data ().end (), value_type/zero/());
400	}
401
402	// Assignment
403	BOOST_UBLAS_INLINE
404	hermitian_matrix &operator = (const hermitian_matrix &m) {
405	size_ = m.size_;
406	data () = m.data ();
407	return *this;
408	}
409	BOOST_UBLAS_INLINE
410	hermitian_matrix &assign_temporary (hermitian_matrix &m) {
411	swap (m);
412	return *this;
413	}
414	template<class AE>
415	BOOST_UBLAS_INLINE
416	hermitian_matrix &operator = (const matrix_expression<AE> &ae) {
417	self_type temporary (ae);
418	return assign_temporary (temporary);
419	}
420	template<class AE>
421	BOOST_UBLAS_INLINE
422	hermitian_matrix &assign (const matrix_expression<AE> &ae) {
423	matrix_assign<scalar_assign> (*this, ae);
424	return *this;
425	}
426	template<class AE>
427	BOOST_UBLAS_INLINE
428	hermitian_matrix& operator += (const matrix_expression<AE> &ae) {
429	self_type temporary (*this + ae);
430	return assign_temporary (temporary);
431	}
432	template<class AE>
433	BOOST_UBLAS_INLINE
434	hermitian_matrix &plus_assign (const matrix_expression<AE> &ae) {
435	matrix_assign<scalar_plus_assign> (*this, ae);
436	return *this;
437	}
438	template<class AE>
439	BOOST_UBLAS_INLINE
440	hermitian_matrix& operator -= (const matrix_expression<AE> &ae) {
441	self_type temporary (*this - ae);
442	return assign_temporary (temporary);
443	}
444	template<class AE>
445	BOOST_UBLAS_INLINE
446	hermitian_matrix &minus_assign (const matrix_expression<AE> &ae) {
447	matrix_assign<scalar_minus_assign> (*this, ae);
448	return *this;
449	}
450	template<class AT>
451	BOOST_UBLAS_INLINE
452	hermitian_matrix& operator *= (const AT &at) {
453	// Multiplication is only allowed for real scalars,
454	// otherwise the resulting matrix isn't hermitian.
455	// Thanks to Peter Schmitteckert for spotting this.
456	BOOST_UBLAS_CHECK (type_traits<value_type>::imag (at) == 0, non_real ());
457	matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
458	return *this;
459	}
460	template<class AT>
461	BOOST_UBLAS_INLINE
462	hermitian_matrix& operator /= (const AT &at) {
463	// Multiplication is only allowed for real scalars,
464	// otherwise the resulting matrix isn't hermitian.
465	// Thanks to Peter Schmitteckert for spotting this.
466	BOOST_UBLAS_CHECK (type_traits<value_type>::imag (at) == 0, non_real ());
467	matrix_assign_scalar<scalar_divides_assign> (*this, at);
468	return *this;
469	}
470
471	// Swapping
472	BOOST_UBLAS_INLINE
473	void swap (hermitian_matrix &m) {
474	if (this != &m) {
475	std::swap (size_, m.size_);
476	data ().swap (m.data ());
477	}
478	}
479	BOOST_UBLAS_INLINE
480	friend void swap (hermitian_matrix &m1, hermitian_matrix &m2) {
481	m1.swap (m2);
482	}
483
484	// Iterator types
485	#ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
486	typedef indexed_iterator1<self_type, packed_random_access_iterator_tag> iterator1;
487	typedef indexed_iterator2<self_type, packed_random_access_iterator_tag> iterator2;
488	typedef indexed_const_iterator1<self_type, packed_random_access_iterator_tag> const_iterator1;
489	typedef indexed_const_iterator2<self_type, packed_random_access_iterator_tag> const_iterator2;
490	#else
491	class const_iterator1;
492	class iterator1;
493	class const_iterator2;
494	class iterator2;
495	#endif
496	typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
497	typedef reverse_iterator_base1<iterator1> reverse_iterator1;
498	typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
499	typedef reverse_iterator_base2<iterator2> reverse_iterator2;
500
501	// Element lookup
502	BOOST_UBLAS_INLINE
503	const_iterator1 find1 (int /* rank */, size_type i, size_type j) const {
504	return const_iterator1 (*this, i, j);
505	}
506	BOOST_UBLAS_INLINE
507	iterator1 find1 (int rank, size_type i, size_type j) {
508	if (rank == 1)
509	i = triangular_type::mutable_restrict1 (i, j);
510	return iterator1 (*this, i, j);
511	}
512	BOOST_UBLAS_INLINE
513	const_iterator2 find2 (int /* rank */, size_type i, size_type j) const {
514	return const_iterator2 (*this, i, j);
515	}
516	BOOST_UBLAS_INLINE
517	iterator2 find2 (int rank, size_type i, size_type j) {
518	if (rank == 1)
519	j = triangular_type::mutable_restrict2 (i, j);
520	return iterator2 (*this, i, j);
521	}
522
523	// Iterators simply are indices.
524
525	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
526	class const_iterator1:
527	public container_const_reference<hermitian_matrix>,
528	public random_access_iterator_base<packed_random_access_iterator_tag,
529	const_iterator1, value_type> {
530	public:
531	typedef typename hermitian_matrix::value_type value_type;
532	typedef typename hermitian_matrix::difference_type difference_type;
533	typedef typename hermitian_matrix::const_reference reference;
534	typedef const typename hermitian_matrix::pointer pointer;
535
536	typedef const_iterator2 dual_iterator_type;
537	typedef const_reverse_iterator2 dual_reverse_iterator_type;
538
539	// Construction and destruction
540	BOOST_UBLAS_INLINE
541	const_iterator1 ():
542	container_const_reference<self_type> (), it1_ (), it2_ () {}
543	BOOST_UBLAS_INLINE
544	const_iterator1 (const self_type &m, size_type it1, size_type it2):
545	container_const_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
546	BOOST_UBLAS_INLINE
547	const_iterator1 (const iterator1 &it):
548	container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
549
550	// Arithmetic
551	BOOST_UBLAS_INLINE
552	const_iterator1 &operator ++ () {
553	++ it1_;
554	return *this;
555	}
556	BOOST_UBLAS_INLINE
557	const_iterator1 &operator -- () {
558	-- it1_;
559	return *this;
560	}
561	BOOST_UBLAS_INLINE
562	const_iterator1 &operator += (difference_type n) {
563	it1_ += n;
564	return *this;
565	}
566	BOOST_UBLAS_INLINE
567	const_iterator1 &operator -= (difference_type n) {
568	it1_ -= n;
569	return *this;
570	}
571	BOOST_UBLAS_INLINE
572	difference_type operator - (const const_iterator1 &it) const {
573	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
574	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
575	return it1_ - it.it1_;
576	}
577
578	// Dereference
579	BOOST_UBLAS_INLINE
580	const_reference operator * () const {
581	return (*this) () (it1_, it2_);
582	}
583
584	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
585	BOOST_UBLAS_INLINE
586	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
587	typename self_type::
588	#endif
589	const_iterator2 begin () const {
590	return (*this) ().find2 (1, it1_, 0);
591	}
592	BOOST_UBLAS_INLINE
593	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
594	typename self_type::
595	#endif
596	const_iterator2 end () const {
597	return (this) ().find2 (1, it1_, (this) ().size2 ());
598	}
599	BOOST_UBLAS_INLINE
600	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
601	typename self_type::
602	#endif
603	const_reverse_iterator2 rbegin () const {
604	return const_reverse_iterator2 (end ());
605	}
606	BOOST_UBLAS_INLINE
607	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
608	typename self_type::
609	#endif
610	const_reverse_iterator2 rend () const {
611	return const_reverse_iterator2 (begin ());
612	}
613	#endif
614
615	// Indices
616	BOOST_UBLAS_INLINE
617	size_type index1 () const {
618	return it1_;
619	}
620	BOOST_UBLAS_INLINE
621	size_type index2 () const {
622	return it2_;
623	}
624
625	// Assignment
626	BOOST_UBLAS_INLINE
627	const_iterator1 &operator = (const const_iterator1 &it) {
628	container_const_reference<self_type>::assign (&it ());
629	it1_ = it.it1_;
630	it2_ = it.it2_;
631	return *this;
632	}
633
634	// Comparison
635	BOOST_UBLAS_INLINE
636	bool operator == (const const_iterator1 &it) const {
637	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
638	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
639	return it1_ == it.it1_;
640	}
641	BOOST_UBLAS_INLINE
642	bool operator < (const const_iterator1 &it) const {
643	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
644	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
645	return it1_ < it.it1_;
646	}
647
648	private:
649	size_type it1_;
650	size_type it2_;
651	};
652	#endif
653
654	BOOST_UBLAS_INLINE
655	const_iterator1 begin1 () const {
656	return find1 (0, 0, 0);
657	}
658	BOOST_UBLAS_INLINE
659	const_iterator1 end1 () const {
660	return find1 (0, size_, 0);
661	}
662
663	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
664	class iterator1:
665	public container_reference<hermitian_matrix>,
666	public random_access_iterator_base<packed_random_access_iterator_tag,
667	iterator1, value_type> {
668	public:
669	typedef typename hermitian_matrix::value_type value_type;
670	typedef typename hermitian_matrix::difference_type difference_type;
671	typedef typename hermitian_matrix::true_reference reference;
672	typedef typename hermitian_matrix::pointer pointer;
673
674	typedef iterator2 dual_iterator_type;
675	typedef reverse_iterator2 dual_reverse_iterator_type;
676
677	// Construction and destruction
678	BOOST_UBLAS_INLINE
679	iterator1 ():
680	container_reference<self_type> (), it1_ (), it2_ () {}
681	BOOST_UBLAS_INLINE
682	iterator1 (self_type &m, size_type it1, size_type it2):
683	container_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
684
685	// Arithmetic
686	BOOST_UBLAS_INLINE
687	iterator1 &operator ++ () {
688	++ it1_;
689	return *this;
690	}
691	BOOST_UBLAS_INLINE
692	iterator1 &operator -- () {
693	-- it1_;
694	return *this;
695	}
696	BOOST_UBLAS_INLINE
697	iterator1 &operator += (difference_type n) {
698	it1_ += n;
699	return *this;
700	}
701	BOOST_UBLAS_INLINE
702	iterator1 &operator -= (difference_type n) {
703	it1_ -= n;
704	return *this;
705	}
706	BOOST_UBLAS_INLINE
707	difference_type operator - (const iterator1 &it) const {
708	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
709	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
710	return it1_ - it.it1_;
711	}
712
713	// Dereference
714	BOOST_UBLAS_INLINE
715	reference operator * () const {
716	return (*this) ().at_element (it1_, it2_);
717	}
718
719	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
720	BOOST_UBLAS_INLINE
721	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
722	typename self_type::
723	#endif
724	iterator2 begin () const {
725	return (*this) ().find2 (1, it1_, 0);
726	}
727	BOOST_UBLAS_INLINE
728	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
729	typename self_type::
730	#endif
731	iterator2 end () const {
732	return (this) ().find2 (1, it1_, (this) ().size2 ());
733	}
734	BOOST_UBLAS_INLINE
735	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
736	typename self_type::
737	#endif
738	reverse_iterator2 rbegin () const {
739	return reverse_iterator2 (end ());
740	}
741	BOOST_UBLAS_INLINE
742	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
743	typename self_type::
744	#endif
745	reverse_iterator2 rend () const {
746	return reverse_iterator2 (begin ());
747	}
748	#endif
749
750	// Indices
751	BOOST_UBLAS_INLINE
752	size_type index1 () const {
753	return it1_;
754	}
755	BOOST_UBLAS_INLINE
756	size_type index2 () const {
757	return it2_;
758	}
759
760	// Assignment
761	BOOST_UBLAS_INLINE
762	iterator1 &operator = (const iterator1 &it) {
763	container_reference<self_type>::assign (&it ());
764	it1_ = it.it1_;
765	it2_ = it.it2_;
766	return *this;
767	}
768
769	// Comparison
770	BOOST_UBLAS_INLINE
771	bool operator == (const iterator1 &it) const {
772	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
773	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
774	return it1_ == it.it1_;
775	}
776	BOOST_UBLAS_INLINE
777	bool operator < (const iterator1 &it) const {
778	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
779	BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
780	return it1_ < it.it1_;
781	}
782
783	private:
784	size_type it1_;
785	size_type it2_;
786
787	friend class const_iterator1;
788	};
789	#endif
790
791	BOOST_UBLAS_INLINE
792	iterator1 begin1 () {
793	return find1 (0, 0, 0);
794	}
795	BOOST_UBLAS_INLINE
796	iterator1 end1 () {
797	return find1 (0, size_, 0);
798	}
799
800	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
801	class const_iterator2:
802	public container_const_reference<hermitian_matrix>,
803	public random_access_iterator_base<packed_random_access_iterator_tag,
804	const_iterator2, value_type> {
805	public:
806	typedef typename hermitian_matrix::value_type value_type;
807	typedef typename hermitian_matrix::difference_type difference_type;
808	typedef typename hermitian_matrix::const_reference reference;
809	typedef const typename hermitian_matrix::pointer pointer;
810
811	typedef const_iterator1 dual_iterator_type;
812	typedef const_reverse_iterator1 dual_reverse_iterator_type;
813
814	// Construction and destruction
815	BOOST_UBLAS_INLINE
816	const_iterator2 ():
817	container_const_reference<self_type> (), it1_ (), it2_ () {}
818	BOOST_UBLAS_INLINE
819	const_iterator2 (const self_type &m, size_type it1, size_type it2):
820	container_const_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
821	BOOST_UBLAS_INLINE
822	const_iterator2 (const iterator2 &it):
823	container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
824
825	// Arithmetic
826	BOOST_UBLAS_INLINE
827	const_iterator2 &operator ++ () {
828	++ it2_;
829	return *this;
830	}
831	BOOST_UBLAS_INLINE
832	const_iterator2 &operator -- () {
833	-- it2_;
834	return *this;
835	}
836	BOOST_UBLAS_INLINE
837	const_iterator2 &operator += (difference_type n) {
838	it2_ += n;
839	return *this;
840	}
841	BOOST_UBLAS_INLINE
842	const_iterator2 &operator -= (difference_type n) {
843	it2_ -= n;
844	return *this;
845	}
846	BOOST_UBLAS_INLINE
847	difference_type operator - (const const_iterator2 &it) const {
848	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
849	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
850	return it2_ - it.it2_;
851	}
852
853	// Dereference
854	BOOST_UBLAS_INLINE
855	const_reference operator * () const {
856	return (*this) () (it1_, it2_);
857	}
858
859	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
860	BOOST_UBLAS_INLINE
861	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
862	typename self_type::
863	#endif
864	const_iterator1 begin () const {
865	return (*this) ().find1 (1, 0, it2_);
866	}
867	BOOST_UBLAS_INLINE
868	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
869	typename self_type::
870	#endif
871	const_iterator1 end () const {
872	return (this) ().find1 (1, (this) ().size1 (), it2_);
873	}
874	BOOST_UBLAS_INLINE
875	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
876	typename self_type::
877	#endif
878	const_reverse_iterator1 rbegin () const {
879	return const_reverse_iterator1 (end ());
880	}
881	BOOST_UBLAS_INLINE
882	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
883	typename self_type::
884	#endif
885	const_reverse_iterator1 rend () const {
886	return const_reverse_iterator1 (begin ());
887	}
888	#endif
889
890	// Indices
891	BOOST_UBLAS_INLINE
892	size_type index1 () const {
893	return it1_;
894	}
895	BOOST_UBLAS_INLINE
896	size_type index2 () const {
897	return it2_;
898	}
899
900	// Assignment
901	BOOST_UBLAS_INLINE
902	const_iterator2 &operator = (const const_iterator2 &it) {
903	container_const_reference<self_type>::assign (&it ());
904	it1_ = it.it1_;
905	it2_ = it.it2_;
906	return *this;
907	}
908
909	// Comparison
910	BOOST_UBLAS_INLINE
911	bool operator == (const const_iterator2 &it) const {
912	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
913	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
914	return it2_ == it.it2_;
915	}
916	BOOST_UBLAS_INLINE
917	bool operator < (const const_iterator2 &it) const {
918	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
919	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
920	return it2_ < it.it2_;
921	}
922
923	private:
924	size_type it1_;
925	size_type it2_;
926	};
927	#endif
928
929	BOOST_UBLAS_INLINE
930	const_iterator2 begin2 () const {
931	return find2 (0, 0, 0);
932	}
933	BOOST_UBLAS_INLINE
934	const_iterator2 end2 () const {
935	return find2 (0, 0, size_);
936	}
937
938	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
939	class iterator2:
940	public container_reference<hermitian_matrix>,
941	public random_access_iterator_base<packed_random_access_iterator_tag,
942	iterator2, value_type> {
943	public:
944	typedef typename hermitian_matrix::value_type value_type;
945	typedef typename hermitian_matrix::difference_type difference_type;
946	typedef typename hermitian_matrix::true_reference reference;
947	typedef typename hermitian_matrix::pointer pointer;
948
949	typedef iterator1 dual_iterator_type;
950	typedef reverse_iterator1 dual_reverse_iterator_type;
951
952	// Construction and destruction
953	BOOST_UBLAS_INLINE
954	iterator2 ():
955	container_reference<self_type> (), it1_ (), it2_ () {}
956	BOOST_UBLAS_INLINE
957	iterator2 (self_type &m, size_type it1, size_type it2):
958	container_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
959
960	// Arithmetic
961	BOOST_UBLAS_INLINE
962	iterator2 &operator ++ () {
963	++ it2_;
964	return *this;
965	}
966	BOOST_UBLAS_INLINE
967	iterator2 &operator -- () {
968	-- it2_;
969	return *this;
970	}
971	BOOST_UBLAS_INLINE
972	iterator2 &operator += (difference_type n) {
973	it2_ += n;
974	return *this;
975	}
976	BOOST_UBLAS_INLINE
977	iterator2 &operator -= (difference_type n) {
978	it2_ -= n;
979	return *this;
980	}
981	BOOST_UBLAS_INLINE
982	difference_type operator - (const iterator2 &it) const {
983	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
984	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
985	return it2_ - it.it2_;
986	}
987
988	// Dereference
989	BOOST_UBLAS_INLINE
990	reference operator * () const {
991	return (*this) ().at_element (it1_, it2_);
992	}
993
994	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
995	BOOST_UBLAS_INLINE
996	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
997	typename self_type::
998	#endif
999	iterator1 begin () const {
1000	return (*this) ().find1 (1, 0, it2_);
1001	}
1002	BOOST_UBLAS_INLINE
1003	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1004	typename self_type::
1005	#endif
1006	iterator1 end () const {
1007	return (this) ().find1 (1, (this) ().size1 (), it2_);
1008	}
1009	BOOST_UBLAS_INLINE
1010	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1011	typename self_type::
1012	#endif
1013	reverse_iterator1 rbegin () const {
1014	return reverse_iterator1 (end ());
1015	}
1016	BOOST_UBLAS_INLINE
1017	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1018	typename self_type::
1019	#endif
1020	reverse_iterator1 rend () const {
1021	return reverse_iterator1 (begin ());
1022	}
1023	#endif
1024
1025	// Indices
1026	BOOST_UBLAS_INLINE
1027	size_type index1 () const {
1028	return it1_;
1029	}
1030	BOOST_UBLAS_INLINE
1031	size_type index2 () const {
1032	return it2_;
1033	}
1034
1035	// Assignment
1036	BOOST_UBLAS_INLINE
1037	iterator2 &operator = (const iterator2 &it) {
1038	container_reference<self_type>::assign (&it ());
1039	it1_ = it.it1_;
1040	it2_ = it.it2_;
1041	return *this;
1042	}
1043
1044	// Comparison
1045	BOOST_UBLAS_INLINE
1046	bool operator == (const iterator2 &it) const {
1047	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1048	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1049	return it2_ == it.it2_;
1050	}
1051	BOOST_UBLAS_INLINE
1052	bool operator < (const iterator2 &it) const {
1053	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1054	BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1055	return it2_ < it.it2_;
1056	}
1057
1058	private:
1059	size_type it1_;
1060	size_type it2_;
1061
1062	friend class const_iterator2;
1063	};
1064	#endif
1065
1066	BOOST_UBLAS_INLINE
1067	iterator2 begin2 () {
1068	return find2 (0, 0, 0);
1069	}
1070	BOOST_UBLAS_INLINE
1071	iterator2 end2 () {
1072	return find2 (0, 0, size_);
1073	}
1074
1075	// Reverse iterators
1076
1077	BOOST_UBLAS_INLINE
1078	const_reverse_iterator1 rbegin1 () const {
1079	return const_reverse_iterator1 (end1 ());
1080	}
1081	BOOST_UBLAS_INLINE
1082	const_reverse_iterator1 rend1 () const {
1083	return const_reverse_iterator1 (begin1 ());
1084	}
1085
1086	BOOST_UBLAS_INLINE
1087	reverse_iterator1 rbegin1 () {
1088	return reverse_iterator1 (end1 ());
1089	}
1090	BOOST_UBLAS_INLINE
1091	reverse_iterator1 rend1 () {
1092	return reverse_iterator1 (begin1 ());
1093	}
1094
1095	BOOST_UBLAS_INLINE
1096	const_reverse_iterator2 rbegin2 () const {
1097	return const_reverse_iterator2 (end2 ());
1098	}
1099	BOOST_UBLAS_INLINE
1100	const_reverse_iterator2 rend2 () const {
1101	return const_reverse_iterator2 (begin2 ());
1102	}
1103
1104	BOOST_UBLAS_INLINE
1105	reverse_iterator2 rbegin2 () {
1106	return reverse_iterator2 (end2 ());
1107	}
1108	BOOST_UBLAS_INLINE
1109	reverse_iterator2 rend2 () {
1110	return reverse_iterator2 (begin2 ());
1111	}
1112
1113	private:
1114	size_type size_;
1115	array_type data_;
1116	};
1117
1118
1119	// Hermitian matrix adaptor class
1120	template<class M, class TRI>
1121	class hermitian_adaptor:
1122	public matrix_expression<hermitian_adaptor<M, TRI> > {
1123
1124	typedef hermitian_adaptor<M, TRI> self_type;
1125	typedef typename M::value_type &true_reference;
1126	public:
1127	#ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1128	using matrix_expression<self_type>::operator ();
1129	#endif
1130	typedef const M const_matrix_type;
1131	typedef M matrix_type;
1132	typedef TRI triangular_type;
1133	typedef typename M::size_type size_type;
1134	typedef typename M::difference_type difference_type;
1135	typedef typename M::value_type value_type;
1136	typedef typename M::value_type const_reference;
1137	#ifndef BOOST_UBLAS_STRICT_HERMITIAN
1138	typedef typename boost::mpl::if_<boost::is_const<M>,
1139	typename M::value_type,
1140	typename M::reference>::type reference;
1141	#else
1142	typedef typename boost::mpl::if_<boost::is_const<M>,
1143	typename M::value_type,
1144	hermitian_matrix_element<self_type> >::type reference;
1145	#endif
1146	typedef typename boost::mpl::if_<boost::is_const<M>,
1147	typename M::const_closure_type,
1148	typename M::closure_type>::type matrix_closure_type;
1149	typedef const self_type const_closure_type;
1150	typedef self_type closure_type;
1151	// Replaced by _temporary_traits to avoid type requirements on M
1152	//typedef typename M::vector_temporary_type vector_temporary_type;
1153	//typedef typename M::matrix_temporary_type matrix_temporary_type;
1154	typedef typename storage_restrict_traits<typename M::storage_category,
1155	packed_proxy_tag>::storage_category storage_category;
1156	typedef typename M::orientation_category orientation_category;
1157
1158	// Construction and destruction
1159	BOOST_UBLAS_INLINE
1160	hermitian_adaptor (matrix_type &data):
1161	matrix_expression<self_type> (),
1162	data_ (data) {
1163	BOOST_UBLAS_CHECK (data_.size1 () == data_.size2 (), bad_size ());
1164	}
1165	BOOST_UBLAS_INLINE
1166	hermitian_adaptor (const hermitian_adaptor &m):
1167	matrix_expression<self_type> (),
1168	data_ (m.data_) {
1169	BOOST_UBLAS_CHECK (data_.size1 () == data_.size2 (), bad_size ());
1170	}
1171
1172	// Accessors
1173	BOOST_UBLAS_INLINE
1174	size_type size1 () const {
1175	return data_.size1 ();
1176	}
1177	BOOST_UBLAS_INLINE
1178	size_type size2 () const {
1179	return data_.size2 ();
1180	}
1181
1182	// Storage accessors
1183	BOOST_UBLAS_INLINE
1184	const matrix_closure_type &data () const {
1185	return data_;
1186	}
1187	BOOST_UBLAS_INLINE
1188	matrix_closure_type &data () {
1189	return data_;
1190	}
1191
1192	// Element access
1193	#ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
1194	BOOST_UBLAS_INLINE
1195	const_reference operator () (size_type i, size_type j) const {
1196	BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1197	BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1198	// if (i == j)
1199	// return type_traits<value_type>::real (data () (i, i));
1200	// else
1201	if (triangular_type::other (i, j))
1202	return data () (i, j);
1203	else
1204	return type_traits<value_type>::conj (data () (j, i));
1205	}
1206	BOOST_UBLAS_INLINE
1207	reference operator () (size_type i, size_type j) {
1208	BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1209	BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1210	#ifndef BOOST_UBLAS_STRICT_HERMITIAN
1211	if (triangular_type::other (i, j))
1212	return data () (i, j);
1213	else {
1214	external_logic ().raise ();
1215	return conj_ = type_traits<value_type>::conj (data () (j, i));
1216	}
1217	#else
1218	if (triangular_type::other (i, j))
1219	return reference (*this, i, j, data () (i, j));
1220	else
1221	return reference (*this, i, j, type_traits<value_type>::conj (data () (j, i)));
1222	#endif
1223	}
1224	BOOST_UBLAS_INLINE
1225	true_reference insert_element (size_type i, size_type j, value_type t) {
1226	BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1227	BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1228	// if (i == j)
1229	// data () (i, i) = type_traits<value_type>::real (t);
1230	// else
1231	if (triangular_type::other (i, j))
1232	return data () (i, j) = t;
1233	else
1234	return data () (j, i) = type_traits<value_type>::conj (t);
1235	}
1236	#else
1237	BOOST_UBLAS_INLINE
1238	reference operator () (size_type i, size_type j) {
1239	BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1240	BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1241	#ifndef BOOST_UBLAS_STRICT_HERMITIAN
1242	if (triangular_type::other (i, j))
1243	return data () (i, j);
1244	else {
1245	external_logic ().raise ();
1246	return conj_ = type_traits<value_type>::conj (data () (j, i));
1247	}
1248	#else
1249	if (triangular_type::other (i, j))
1250	return reference (*this, i, j, data () (i, j));
1251	else
1252	return reference (*this, i, j, type_traits<value_type>::conj (data () (j, i)));
1253	#endif
1254	}
1255	BOOST_UBLAS_INLINE
1256	true_reference insert_element (size_type i, size_type j, value_type t) {
1257	BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1258	BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1259	// if (i == j)
1260	// data () (i, i) = type_traits<value_type>::real (t);
1261	// else
1262	if (triangular_type::other (i, j))
1263	return data () (i, j) = t;
1264	else
1265	return data () (j, i) = type_traits<value_type>::conj (t);
1266	}
1267	#endif
1268
1269	// Assignment
1270	BOOST_UBLAS_INLINE
1271	hermitian_adaptor &operator = (const hermitian_adaptor &m) {
1272	matrix_assign<scalar_assign, triangular_type> (*this, m);
1273	return *this;
1274	}
1275	BOOST_UBLAS_INLINE
1276	hermitian_adaptor &assign_temporary (hermitian_adaptor &m) {
1277	*this = m;
1278	return *this;
1279	}
1280	template<class AE>
1281	BOOST_UBLAS_INLINE
1282	hermitian_adaptor &operator = (const matrix_expression<AE> &ae) {
1283	matrix_assign<scalar_assign, triangular_type> (*this, matrix<value_type> (ae));
1284	return *this;
1285	}
1286	template<class AE>
1287	BOOST_UBLAS_INLINE
1288	hermitian_adaptor &assign (const matrix_expression<AE> &ae) {
1289	matrix_assign<scalar_assign, triangular_type> (*this, ae);
1290	return *this;
1291	}
1292	template<class AE>
1293	BOOST_UBLAS_INLINE
1294	hermitian_adaptor& operator += (const matrix_expression<AE> &ae) {
1295	matrix_assign<scalar_assign, triangular_type> (this, matrix<value_type> (this + ae));
1296	return *this;
1297	}
1298	template<class AE>
1299	BOOST_UBLAS_INLINE
1300	hermitian_adaptor &plus_assign (const matrix_expression<AE> &ae) {
1301	matrix_assign<scalar_plus_assign, triangular_type> (*this, ae);
1302	return *this;
1303	}
1304	template<class AE>
1305	BOOST_UBLAS_INLINE
1306	hermitian_adaptor& operator -= (const matrix_expression<AE> &ae) {
1307	matrix_assign<scalar_assign, triangular_type> (this, matrix<value_type> (this - ae));
1308	return *this;
1309	}
1310	template<class AE>
1311	BOOST_UBLAS_INLINE
1312	hermitian_adaptor &minus_assign (const matrix_expression<AE> &ae) {
1313	matrix_assign<scalar_minus_assign, triangular_type> (*this, ae);
1314	return *this;
1315	}
1316	template<class AT>
1317	BOOST_UBLAS_INLINE
1318	hermitian_adaptor& operator *= (const AT &at) {
1319	// Multiplication is only allowed for real scalars,
1320	// otherwise the resulting matrix isn't hermitian.
1321	// Thanks to Peter Schmitteckert for spotting this.
1322	BOOST_UBLAS_CHECK (type_traits<value_type>::imag (at) == 0, non_real ());
1323	matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
1324	return *this;
1325	}
1326	template<class AT>
1327	BOOST_UBLAS_INLINE
1328	hermitian_adaptor& operator /= (const AT &at) {
1329	// Multiplication is only allowed for real scalars,
1330	// otherwise the resulting matrix isn't hermitian.
1331	// Thanks to Peter Schmitteckert for spotting this.
1332	BOOST_UBLAS_CHECK (type_traits<value_type>::imag (at) == 0, non_real ());
1333	matrix_assign_scalar<scalar_divides_assign> (*this, at);
1334	return *this;
1335	}
1336
1337	// Closure comparison
1338	BOOST_UBLAS_INLINE
1339	bool same_closure (const hermitian_adaptor &ha) const {
1340	return (*this).data ().same_closure (ha.data ());
1341	}
1342
1343	// Swapping
1344	BOOST_UBLAS_INLINE
1345	void swap (hermitian_adaptor &m) {
1346	if (this != &m)
1347	matrix_swap<scalar_swap, triangular_type> (*this, m);
1348	}
1349	BOOST_UBLAS_INLINE
1350	friend void swap (hermitian_adaptor &m1, hermitian_adaptor &m2) {
1351	m1.swap (m2);
1352	}
1353
1354	// Iterator types
1355	private:
1356	// Use matrix iterator
1357	typedef typename M::const_iterator1 const_subiterator1_type;
1358	typedef typename boost::mpl::if_<boost::is_const<M>,
1359	typename M::const_iterator1,
1360	typename M::iterator1>::type subiterator1_type;
1361	typedef typename M::const_iterator2 const_subiterator2_type;
1362	typedef typename boost::mpl::if_<boost::is_const<M>,
1363	typename M::const_iterator2,
1364	typename M::iterator2>::type subiterator2_type;
1365
1366	public:
1367	#ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1368	typedef indexed_iterator1<self_type, packed_random_access_iterator_tag> iterator1;
1369	typedef indexed_iterator2<self_type, packed_random_access_iterator_tag> iterator2;
1370	typedef indexed_const_iterator1<self_type, dense_random_access_iterator_tag> const_iterator1;
1371	typedef indexed_const_iterator2<self_type, dense_random_access_iterator_tag> const_iterator2;
1372	#else
1373	class const_iterator1;
1374	class iterator1;
1375	class const_iterator2;
1376	class iterator2;
1377	#endif
1378	typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
1379	typedef reverse_iterator_base1<iterator1> reverse_iterator1;
1380	typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
1381	typedef reverse_iterator_base2<iterator2> reverse_iterator2;
1382
1383	// Element lookup
1384	BOOST_UBLAS_INLINE
1385	const_iterator1 find1 (int rank, size_type i, size_type j) const {
1386	if (triangular_type::other (i, j)) {
1387	if (triangular_type::other (size1 (), j)) {
1388	return const_iterator1 (*this, 0, 0,
1389	data ().find1 (rank, i, j), data ().find1 (rank, size1 (), j),
1390	data ().find2 (rank, size2 (), size1 ()), data ().find2 (rank, size2 (), size1 ()));
1391	} else {
1392	return const_iterator1 (*this, 0, 1,
1393	data ().find1 (rank, i, j), data ().find1 (rank, j, j),
1394	data ().find2 (rank, j, j), data ().find2 (rank, j, size1 ()));
1395	}
1396	} else {
1397	if (triangular_type::other (size1 (), j)) {
1398	return const_iterator1 (*this, 1, 0,
1399	data ().find1 (rank, j, j), data ().find1 (rank, size1 (), j),
1400	data ().find2 (rank, j, i), data ().find2 (rank, j, j));
1401	} else {
1402	return const_iterator1 (*this, 1, 1,
1403	data ().find1 (rank, size1 (), size2 ()), data ().find1 (rank, size1 (), size2 ()),
1404	data ().find2 (rank, j, i), data ().find2 (rank, j, size1 ()));
1405	}
1406	}
1407	}
1408	BOOST_UBLAS_INLINE
1409	iterator1 find1 (int rank, size_type i, size_type j) {
1410	if (rank == 1)
1411	i = triangular_type::mutable_restrict1 (i, j);
1412	return iterator1 (*this, data ().find1 (rank, i, j));
1413	}
1414	BOOST_UBLAS_INLINE
1415	const_iterator2 find2 (int rank, size_type i, size_type j) const {
1416	if (triangular_type::other (i, j)) {
1417	if (triangular_type::other (i, size2 ())) {
1418	return const_iterator2 (*this, 1, 1,
1419	data ().find1 (rank, size2 (), size1 ()), data ().find1 (rank, size2 (), size1 ()),
1420	data ().find2 (rank, i, j), data ().find2 (rank, i, size2 ()));
1421	} else {
1422	return const_iterator2 (*this, 1, 0,
1423	data ().find1 (rank, i, i), data ().find1 (rank, size2 (), i),
1424	data ().find2 (rank, i, j), data ().find2 (rank, i, i));
1425	}
1426	} else {
1427	if (triangular_type::other (i, size2 ())) {
1428	return const_iterator2 (*this, 0, 1,
1429	data ().find1 (rank, j, i), data ().find1 (rank, i, i),
1430	data ().find2 (rank, i, i), data ().find2 (rank, i, size2 ()));
1431	} else {
1432	return const_iterator2 (*this, 0, 0,
1433	data ().find1 (rank, j, i), data ().find1 (rank, size2 (), i),
1434	data ().find2 (rank, size1 (), size2 ()), data ().find2 (rank, size2 (), size2 ()));
1435	}
1436	}
1437	}
1438	BOOST_UBLAS_INLINE
1439	iterator2 find2 (int rank, size_type i, size_type j) {
1440	if (rank == 1)
1441	j = triangular_type::mutable_restrict2 (i, j);
1442	return iterator2 (*this, data ().find2 (rank, i, j));
1443	}
1444
1445	// Iterators simply are indices.
1446
1447	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1448	class const_iterator1:
1449	public container_const_reference<hermitian_adaptor>,
1450	public random_access_iterator_base<typename iterator_restrict_traits<
1451	typename const_subiterator1_type::iterator_category, dense_random_access_iterator_tag>::iterator_category,
1452	const_iterator1, value_type> {
1453	public:
1454	typedef typename const_subiterator1_type::value_type value_type;
1455	typedef typename const_subiterator1_type::difference_type difference_type;
1456	// FIXME no better way to not return the address of a temporary?
1457	// typedef typename const_subiterator1_type::reference reference;
1458	typedef typename const_subiterator1_type::value_type reference;
1459	typedef typename const_subiterator1_type::pointer pointer;
1460
1461	typedef const_iterator2 dual_iterator_type;
1462	typedef const_reverse_iterator2 dual_reverse_iterator_type;
1463
1464	// Construction and destruction
1465	BOOST_UBLAS_INLINE
1466	const_iterator1 ():
1467	container_const_reference<self_type> (),
1468	begin_ (-1), end_ (-1), current_ (-1),
1469	it1_begin_ (), it1_end_ (), it1_ (),
1470	it2_begin_ (), it2_end_ (), it2_ () {}
1471	BOOST_UBLAS_INLINE
1472	const_iterator1 (const self_type &m, int begin, int end,
1473	const const_subiterator1_type &it1_begin, const const_subiterator1_type &it1_end,
1474	const const_subiterator2_type &it2_begin, const const_subiterator2_type &it2_end):
1475	container_const_reference<self_type> (m),
1476	begin_ (begin), end_ (end), current_ (begin),
1477	it1_begin_ (it1_begin), it1_end_ (it1_end), it1_ (it1_begin_),
1478	it2_begin_ (it2_begin), it2_end_ (it2_end), it2_ (it2_begin_) {
1479	if (current_ == 0 && it1_ == it1_end_)
1480	current_ = 1;
1481	if (current_ == 1 && it2_ == it2_end_)
1482	current_ = 0;
1483	if ((current_ == 0 && it1_ == it1_end_) \|\|
1484	(current_ == 1 && it2_ == it2_end_))
1485	current_ = end_;
1486	BOOST_UBLAS_CHECK (current_ == end_ \|\|
1487	(current_ == 0 && it1_ != it1_end_) \|\|
1488	(current_ == 1 && it2_ != it2_end_), internal_logic ());
1489	}
1490	// FIXME cannot compile
1491	// iterator1 does not have these members!
1492	BOOST_UBLAS_INLINE
1493	const_iterator1 (const iterator1 &it):
1494	container_const_reference<self_type> (it ()),
1495	begin_ (it.begin_), end_ (it.end_), current_ (it.current_),
1496	it1_begin_ (it.it1_begin_), it1_end_ (it.it1_end_), it1_ (it.it1_),
1497	it2_begin_ (it.it2_begin_), it2_end_ (it.it2_end_), it2_ (it.it2_) {
1498	BOOST_UBLAS_CHECK (current_ == end_ \|\|
1499	(current_ == 0 && it1_ != it1_end_) \|\|
1500	(current_ == 1 && it2_ != it2_end_), internal_logic ());
1501	}
1502
1503	// Arithmetic
1504	BOOST_UBLAS_INLINE
1505	const_iterator1 &operator ++ () {
1506	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1507	if (current_ == 0) {
1508	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
1509	++ it1_;
1510	if (it1_ == it1_end_ && end_ == 1) {
1511	it2_ = it2_begin_;
1512	current_ = 1;
1513	}
1514	} else /* if (current_ == 1) */ {
1515	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
1516	++ it2_;
1517	if (it2_ == it2_end_ && end_ == 0) {
1518	it1_ = it1_begin_;
1519	current_ = 0;
1520	}
1521	}
1522	return *this;
1523	}
1524	BOOST_UBLAS_INLINE
1525	const_iterator1 &operator -- () {
1526	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1527	if (current_ == 0) {
1528	if (it1_ == it1_begin_ && begin_ == 1) {
1529	it2_ = it2_end_;
1530	BOOST_UBLAS_CHECK (it2_ != it2_begin_, internal_logic ());
1531	-- it2_;
1532	current_ = 1;
1533	} else {
1534	-- it1_;
1535	}
1536	} else /* if (current_ == 1) */ {
1537	if (it2_ == it2_begin_ && begin_ == 0) {
1538	it1_ = it1_end_;
1539	BOOST_UBLAS_CHECK (it1_ != it1_begin_, internal_logic ());
1540	-- it1_;
1541	current_ = 0;
1542	} else {
1543	-- it2_;
1544	}
1545	}
1546	return *this;
1547	}
1548	BOOST_UBLAS_INLINE
1549	const_iterator1 &operator += (difference_type n) {
1550	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1551	if (current_ == 0) {
1552	size_type d = (std::min) (n, it1_end_ - it1_);
1553	it1_ += d;
1554	n -= d;
1555	if (n > 0 \|\| (end_ == 1 && it1_ == it1_end_)) {
1556	BOOST_UBLAS_CHECK (end_ == 1, external_logic ());
1557	d = (std::min) (n, it2_end_ - it2_begin_);
1558	it2_ = it2_begin_ + d;
1559	n -= d;
1560	current_ = 1;
1561	}
1562	} else /* if (current_ == 1) */ {
1563	size_type d = (std::min) (n, it2_end_ - it2_);
1564	it2_ += d;
1565	n -= d;
1566	if (n > 0 \|\| (end_ == 0 && it2_ == it2_end_)) {
1567	BOOST_UBLAS_CHECK (end_ == 0, external_logic ());
1568	d = (std::min) (n, it1_end_ - it1_begin_);
1569	it1_ = it1_begin_ + d;
1570	n -= d;
1571	current_ = 0;
1572	}
1573	}
1574	BOOST_UBLAS_CHECK (n == 0, external_logic ());
1575	return *this;
1576	}
1577	BOOST_UBLAS_INLINE
1578	const_iterator1 &operator -= (difference_type n) {
1579	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1580	if (current_ == 0) {
1581	size_type d = (std::min) (n, it1_ - it1_begin_);
1582	it1_ -= d;
1583	n -= d;
1584	if (n > 0) {
1585	BOOST_UBLAS_CHECK (end_ == 1, external_logic ());
1586	d = (std::min) (n, it2_end_ - it2_begin_);
1587	it2_ = it2_end_ - d;
1588	n -= d;
1589	current_ = 1;
1590	}
1591	} else /* if (current_ == 1) */ {
1592	size_type d = (std::min) (n, it2_ - it2_begin_);
1593	it2_ -= d;
1594	n -= d;
1595	if (n > 0) {
1596	BOOST_UBLAS_CHECK (end_ == 0, external_logic ());
1597	d = (std::min) (n, it1_end_ - it1_begin_);
1598	it1_ = it1_end_ - d;
1599	n -= d;
1600	current_ = 0;
1601	}
1602	}
1603	BOOST_UBLAS_CHECK (n == 0, external_logic ());
1604	return *this;
1605	}
1606	BOOST_UBLAS_INLINE
1607	difference_type operator - (const const_iterator1 &it) const {
1608	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1609	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1610	BOOST_UBLAS_CHECK (it.current_ == 0 \|\| it.current_ == 1, internal_logic ());
1611	BOOST_UBLAS_CHECK (/* begin_ == it.begin_ && */ end_ == it.end_, internal_logic ());
1612	if (current_ == 0 && it.current_ == 0) {
1613	return it1_ - it.it1_;
1614	} else if (current_ == 0 && it.current_ == 1) {
1615	if (end_ == 1 && it.end_ == 1) {
1616	return (it1_ - it.it1_end_) + (it.it2_begin_ - it.it2_);
1617	} else /* if (end_ == 0 && it.end_ == 0) */ {
1618	return (it1_ - it.it1_begin_) + (it.it2_end_ - it.it2_);
1619	}
1620
1621	} else if (current_ == 1 && it.current_ == 0) {
1622	if (end_ == 1 && it.end_ == 1) {
1623	return (it2_ - it.it2_begin_) + (it.it1_end_ - it.it1_);
1624	} else /* if (end_ == 0 && it.end_ == 0) */ {
1625	return (it2_ - it.it2_end_) + (it.it1_begin_ - it.it1_);
1626	}
1627	} else /* if (current_ == 1 && it.current_ == 1) */ {
1628	return it2_ - it.it2_;
1629	}
1630	}
1631
1632	// Dereference
1633	BOOST_UBLAS_INLINE
1634	const_reference operator * () const {
1635	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1636	if (current_ == 0) {
1637	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
1638	if (triangular_type::other (index1 (), index2 ()))
1639	return *it1_;
1640	else
1641	return type_traits<value_type>::conj (*it1_);
1642	} else /* if (current_ == 1) */ {
1643	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
1644	if (triangular_type::other (index1 (), index2 ()))
1645	return *it2_;
1646	else
1647	return type_traits<value_type>::conj (*it2_);
1648	}
1649	}
1650
1651	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1652	BOOST_UBLAS_INLINE
1653	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1654	typename self_type::
1655	#endif
1656	const_iterator2 begin () const {
1657	return (*this) ().find2 (1, index1 (), 0);
1658	}
1659	BOOST_UBLAS_INLINE
1660	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1661	typename self_type::
1662	#endif
1663	const_iterator2 end () const {
1664	return (this) ().find2 (1, index1 (), (this) ().size2 ());
1665	}
1666	BOOST_UBLAS_INLINE
1667	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1668	typename self_type::
1669	#endif
1670	const_reverse_iterator2 rbegin () const {
1671	return const_reverse_iterator2 (end ());
1672	}
1673	BOOST_UBLAS_INLINE
1674	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1675	typename self_type::
1676	#endif
1677	const_reverse_iterator2 rend () const {
1678	return const_reverse_iterator2 (begin ());
1679	}
1680	#endif
1681
1682	// Indices
1683	BOOST_UBLAS_INLINE
1684	size_type index1 () const {
1685	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1686	if (current_ == 0) {
1687	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
1688	return it1_.index1 ();
1689	} else /* if (current_ == 1) */ {
1690	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
1691	return it2_.index2 ();
1692	}
1693	}
1694	BOOST_UBLAS_INLINE
1695	size_type index2 () const {
1696	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1697	if (current_ == 0) {
1698	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
1699	return it1_.index2 ();
1700	} else /* if (current_ == 1) */ {
1701	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
1702	return it2_.index1 ();
1703	}
1704	}
1705
1706	// Assignment
1707	BOOST_UBLAS_INLINE
1708	const_iterator1 &operator = (const const_iterator1 &it) {
1709	container_const_reference<self_type>::assign (&it ());
1710	begin_ = it.begin_;
1711	end_ = it.end_;
1712	current_ = it.current_;
1713	it1_begin_ = it.it1_begin_;
1714	it1_end_ = it.it1_end_;
1715	it1_ = it.it1_;
1716	it2_begin_ = it.it2_begin_;
1717	it2_end_ = it.it2_end_;
1718	it2_ = it.it2_;
1719	return *this;
1720	}
1721
1722	// Comparison
1723	BOOST_UBLAS_INLINE
1724	bool operator == (const const_iterator1 &it) const {
1725	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1726	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1727	BOOST_UBLAS_CHECK (it.current_ == 0 \|\| it.current_ == 1, internal_logic ());
1728	BOOST_UBLAS_CHECK (/* begin_ == it.begin_ && */ end_ == it.end_, internal_logic ());
1729	return (current_ == 0 && it.current_ == 0 && it1_ == it.it1_) \|\|
1730	(current_ == 1 && it.current_ == 1 && it2_ == it.it2_);
1731	}
1732	BOOST_UBLAS_INLINE
1733	bool operator < (const const_iterator1 &it) const {
1734	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1735	return it - *this > 0;
1736	}
1737
1738	private:
1739	int begin_;
1740	int end_;
1741	int current_;
1742	const_subiterator1_type it1_begin_;
1743	const_subiterator1_type it1_end_;
1744	const_subiterator1_type it1_;
1745	const_subiterator2_type it2_begin_;
1746	const_subiterator2_type it2_end_;
1747	const_subiterator2_type it2_;
1748	};
1749	#endif
1750
1751	BOOST_UBLAS_INLINE
1752	const_iterator1 begin1 () const {
1753	return find1 (0, 0, 0);
1754	}
1755	BOOST_UBLAS_INLINE
1756	const_iterator1 end1 () const {
1757	return find1 (0, size1 (), 0);
1758	}
1759
1760	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1761	class iterator1:
1762	public container_reference<hermitian_adaptor>,
1763	public random_access_iterator_base<typename iterator_restrict_traits<
1764	typename subiterator1_type::iterator_category, packed_random_access_iterator_tag>::iterator_category,
1765	iterator1, value_type> {
1766	public:
1767	typedef typename subiterator1_type::value_type value_type;
1768	typedef typename subiterator1_type::difference_type difference_type;
1769	typedef typename subiterator1_type::reference reference;
1770	typedef typename subiterator1_type::pointer pointer;
1771
1772	typedef iterator2 dual_iterator_type;
1773	typedef reverse_iterator2 dual_reverse_iterator_type;
1774
1775	// Construction and destruction
1776	BOOST_UBLAS_INLINE
1777	iterator1 ():
1778	container_reference<self_type> (), it1_ () {}
1779	BOOST_UBLAS_INLINE
1780	iterator1 (self_type &m, const subiterator1_type &it1):
1781	container_reference<self_type> (m), it1_ (it1) {}
1782
1783	// Arithmetic
1784	BOOST_UBLAS_INLINE
1785	iterator1 &operator ++ () {
1786	++ it1_;
1787	return *this;
1788	}
1789	BOOST_UBLAS_INLINE
1790	iterator1 &operator -- () {
1791	-- it1_;
1792	return *this;
1793	}
1794	BOOST_UBLAS_INLINE
1795	iterator1 &operator += (difference_type n) {
1796	it1_ += n;
1797	return *this;
1798	}
1799	BOOST_UBLAS_INLINE
1800	iterator1 &operator -= (difference_type n) {
1801	it1_ -= n;
1802	return *this;
1803	}
1804	BOOST_UBLAS_INLINE
1805	difference_type operator - (const iterator1 &it) const {
1806	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1807	return it1_ - it.it1_;
1808	}
1809
1810	// Dereference
1811	BOOST_UBLAS_INLINE
1812	reference operator * () const {
1813	return *it1_;
1814	}
1815
1816	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1817	BOOST_UBLAS_INLINE
1818	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1819	typename self_type::
1820	#endif
1821	iterator2 begin () const {
1822	return (*this) ().find2 (1, index1 (), 0);
1823	}
1824	BOOST_UBLAS_INLINE
1825	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1826	typename self_type::
1827	#endif
1828	iterator2 end () const {
1829	return (this) ().find2 (1, index1 (), (this) ().size2 ());
1830	}
1831	BOOST_UBLAS_INLINE
1832	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1833	typename self_type::
1834	#endif
1835	reverse_iterator2 rbegin () const {
1836	return reverse_iterator2 (end ());
1837	}
1838	BOOST_UBLAS_INLINE
1839	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1840	typename self_type::
1841	#endif
1842	reverse_iterator2 rend () const {
1843	return reverse_iterator2 (begin ());
1844	}
1845	#endif
1846
1847	// Indices
1848	BOOST_UBLAS_INLINE
1849	size_type index1 () const {
1850	return it1_.index1 ();
1851	}
1852	BOOST_UBLAS_INLINE
1853	size_type index2 () const {
1854	return it1_.index2 ();
1855	}
1856
1857	// Assignment
1858	BOOST_UBLAS_INLINE
1859	iterator1 &operator = (const iterator1 &it) {
1860	container_reference<self_type>::assign (&it ());
1861	it1_ = it.it1_;
1862	return *this;
1863	}
1864
1865	// Comparison
1866	BOOST_UBLAS_INLINE
1867	bool operator == (const iterator1 &it) const {
1868	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1869	return it1_ == it.it1_;
1870	}
1871	BOOST_UBLAS_INLINE
1872	bool operator < (const iterator1 &it) const {
1873	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1874	return it1_ < it.it1_;
1875	}
1876
1877	private:
1878	subiterator1_type it1_;
1879
1880	friend class const_iterator1;
1881	};
1882	#endif
1883
1884	BOOST_UBLAS_INLINE
1885	iterator1 begin1 () {
1886	return find1 (0, 0, 0);
1887	}
1888	BOOST_UBLAS_INLINE
1889	iterator1 end1 () {
1890	return find1 (0, size1 (), 0);
1891	}
1892
1893	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1894	class const_iterator2:
1895	public container_const_reference<hermitian_adaptor>,
1896	public random_access_iterator_base<typename iterator_restrict_traits<
1897	typename const_subiterator2_type::iterator_category, dense_random_access_iterator_tag>::iterator_category,
1898	const_iterator2, value_type> {
1899	public:
1900	typedef typename const_subiterator2_type::value_type value_type;
1901	typedef typename const_subiterator2_type::difference_type difference_type;
1902	// FIXME no better way to not return the address of a temporary?
1903	// typedef typename const_subiterator2_type::reference reference;
1904	typedef typename const_subiterator2_type::value_type reference;
1905	typedef typename const_subiterator2_type::pointer pointer;
1906
1907	typedef const_iterator1 dual_iterator_type;
1908	typedef const_reverse_iterator1 dual_reverse_iterator_type;
1909
1910	// Construction and destruction
1911	BOOST_UBLAS_INLINE
1912	const_iterator2 ():
1913	container_const_reference<self_type> (),
1914	begin_ (-1), end_ (-1), current_ (-1),
1915	it1_begin_ (), it1_end_ (), it1_ (),
1916	it2_begin_ (), it2_end_ (), it2_ () {}
1917	BOOST_UBLAS_INLINE
1918	const_iterator2 (const self_type &m, int begin, int end,
1919	const const_subiterator1_type &it1_begin, const const_subiterator1_type &it1_end,
1920	const const_subiterator2_type &it2_begin, const const_subiterator2_type &it2_end):
1921	container_const_reference<self_type> (m),
1922	begin_ (begin), end_ (end), current_ (begin),
1923	it1_begin_ (it1_begin), it1_end_ (it1_end), it1_ (it1_begin_),
1924	it2_begin_ (it2_begin), it2_end_ (it2_end), it2_ (it2_begin_) {
1925	if (current_ == 0 && it1_ == it1_end_)
1926	current_ = 1;
1927	if (current_ == 1 && it2_ == it2_end_)
1928	current_ = 0;
1929	if ((current_ == 0 && it1_ == it1_end_) \|\|
1930	(current_ == 1 && it2_ == it2_end_))
1931	current_ = end_;
1932	BOOST_UBLAS_CHECK (current_ == end_ \|\|
1933	(current_ == 0 && it1_ != it1_end_) \|\|
1934	(current_ == 1 && it2_ != it2_end_), internal_logic ());
1935	}
1936	// FIXME cannot compiler
1937	// iterator2 does not have these members!
1938	BOOST_UBLAS_INLINE
1939	const_iterator2 (const iterator2 &it):
1940	container_const_reference<self_type> (it ()),
1941	begin_ (it.begin_), end_ (it.end_), current_ (it.current_),
1942	it1_begin_ (it.it1_begin_), it1_end_ (it.it1_end_), it1_ (it.it1_),
1943	it2_begin_ (it.it2_begin_), it2_end_ (it.it2_end_), it2_ (it.it2_) {
1944	BOOST_UBLAS_CHECK (current_ == end_ \|\|
1945	(current_ == 0 && it1_ != it1_end_) \|\|
1946	(current_ == 1 && it2_ != it2_end_), internal_logic ());
1947	}
1948
1949	// Arithmetic
1950	BOOST_UBLAS_INLINE
1951	const_iterator2 &operator ++ () {
1952	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1953	if (current_ == 0) {
1954	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
1955	++ it1_;
1956	if (it1_ == it1_end_ && end_ == 1) {
1957	it2_ = it2_begin_;
1958	current_ = 1;
1959	}
1960	} else /* if (current_ == 1) */ {
1961	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
1962	++ it2_;
1963	if (it2_ == it2_end_ && end_ == 0) {
1964	it1_ = it1_begin_;
1965	current_ = 0;
1966	}
1967	}
1968	return *this;
1969	}
1970	BOOST_UBLAS_INLINE
1971	const_iterator2 &operator -- () {
1972	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1973	if (current_ == 0) {
1974	if (it1_ == it1_begin_ && begin_ == 1) {
1975	it2_ = it2_end_;
1976	BOOST_UBLAS_CHECK (it2_ != it2_begin_, internal_logic ());
1977	-- it2_;
1978	current_ = 1;
1979	} else {
1980	-- it1_;
1981	}
1982	} else /* if (current_ == 1) */ {
1983	if (it2_ == it2_begin_ && begin_ == 0) {
1984	it1_ = it1_end_;
1985	BOOST_UBLAS_CHECK (it1_ != it1_begin_, internal_logic ());
1986	-- it1_;
1987	current_ = 0;
1988	} else {
1989	-- it2_;
1990	}
1991	}
1992	return *this;
1993	}
1994	BOOST_UBLAS_INLINE
1995	const_iterator2 &operator += (difference_type n) {
1996	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
1997	if (current_ == 0) {
1998	size_type d = (std::min) (n, it1_end_ - it1_);
1999	it1_ += d;
2000	n -= d;
2001	if (n > 0 \|\| (end_ == 1 && it1_ == it1_end_)) {
2002	BOOST_UBLAS_CHECK (end_ == 1, external_logic ());
2003	d = (std::min) (n, it2_end_ - it2_begin_);
2004	it2_ = it2_begin_ + d;
2005	n -= d;
2006	current_ = 1;
2007	}
2008	} else /* if (current_ == 1) */ {
2009	size_type d = (std::min) (n, it2_end_ - it2_);
2010	it2_ += d;
2011	n -= d;
2012	if (n > 0 \|\| (end_ == 0 && it2_ == it2_end_)) {
2013	BOOST_UBLAS_CHECK (end_ == 0, external_logic ());
2014	d = (std::min) (n, it1_end_ - it1_begin_);
2015	it1_ = it1_begin_ + d;
2016	n -= d;
2017	current_ = 0;
2018	}
2019	}
2020	BOOST_UBLAS_CHECK (n == 0, external_logic ());
2021	return *this;
2022	}
2023	BOOST_UBLAS_INLINE
2024	const_iterator2 &operator -= (difference_type n) {
2025	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2026	if (current_ == 0) {
2027	size_type d = (std::min) (n, it1_ - it1_begin_);
2028	it1_ -= d;
2029	n -= d;
2030	if (n > 0) {
2031	BOOST_UBLAS_CHECK (end_ == 1, external_logic ());
2032	d = (std::min) (n, it2_end_ - it2_begin_);
2033	it2_ = it2_end_ - d;
2034	n -= d;
2035	current_ = 1;
2036	}
2037	} else /* if (current_ == 1) */ {
2038	size_type d = (std::min) (n, it2_ - it2_begin_);
2039	it2_ -= d;
2040	n -= d;
2041	if (n > 0) {
2042	BOOST_UBLAS_CHECK (end_ == 0, external_logic ());
2043	d = (std::min) (n, it1_end_ - it1_begin_);
2044	it1_ = it1_end_ - d;
2045	n -= d;
2046	current_ = 0;
2047	}
2048	}
2049	BOOST_UBLAS_CHECK (n == 0, external_logic ());
2050	return *this;
2051	}
2052	BOOST_UBLAS_INLINE
2053	difference_type operator - (const const_iterator2 &it) const {
2054	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2055	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2056	BOOST_UBLAS_CHECK (it.current_ == 0 \|\| it.current_ == 1, internal_logic ());
2057	BOOST_UBLAS_CHECK (/* begin_ == it.begin_ && */ end_ == it.end_, internal_logic ());
2058	if (current_ == 0 && it.current_ == 0) {
2059	return it1_ - it.it1_;
2060	} else if (current_ == 0 && it.current_ == 1) {
2061	if (end_ == 1 && it.end_ == 1) {
2062	return (it1_ - it.it1_end_) + (it.it2_begin_ - it.it2_);
2063	} else /* if (end_ == 0 && it.end_ == 0) */ {
2064	return (it1_ - it.it1_begin_) + (it.it2_end_ - it.it2_);
2065	}
2066
2067	} else if (current_ == 1 && it.current_ == 0) {
2068	if (end_ == 1 && it.end_ == 1) {
2069	return (it2_ - it.it2_begin_) + (it.it1_end_ - it.it1_);
2070	} else /* if (end_ == 0 && it.end_ == 0) */ {
2071	return (it2_ - it.it2_end_) + (it.it1_begin_ - it.it1_);
2072	}
2073	} else /* if (current_ == 1 && it.current_ == 1) */ {
2074	return it2_ - it.it2_;
2075	}
2076	}
2077
2078	// Dereference
2079	BOOST_UBLAS_INLINE
2080	const_reference operator * () const {
2081	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2082	if (current_ == 0) {
2083	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
2084	if (triangular_type::other (index1 (), index2 ()))
2085	return *it1_;
2086	else
2087	return type_traits<value_type>::conj (*it1_);
2088	} else /* if (current_ == 1) */ {
2089	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
2090	if (triangular_type::other (index1 (), index2 ()))
2091	return *it2_;
2092	else
2093	return type_traits<value_type>::conj (*it2_);
2094	}
2095	}
2096
2097	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2098	BOOST_UBLAS_INLINE
2099	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2100	typename self_type::
2101	#endif
2102	const_iterator1 begin () const {
2103	return (*this) ().find1 (1, 0, index2 ());
2104	}
2105	BOOST_UBLAS_INLINE
2106	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2107	typename self_type::
2108	#endif
2109	const_iterator1 end () const {
2110	return (this) ().find1 (1, (this) ().size1 (), index2 ());
2111	}
2112	BOOST_UBLAS_INLINE
2113	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2114	typename self_type::
2115	#endif
2116	const_reverse_iterator1 rbegin () const {
2117	return const_reverse_iterator1 (end ());
2118	}
2119	BOOST_UBLAS_INLINE
2120	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2121	typename self_type::
2122	#endif
2123	const_reverse_iterator1 rend () const {
2124	return const_reverse_iterator1 (begin ());
2125	}
2126	#endif
2127
2128	// Indices
2129	BOOST_UBLAS_INLINE
2130	size_type index1 () const {
2131	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2132	if (current_ == 0) {
2133	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
2134	return it1_.index2 ();
2135	} else /* if (current_ == 1) */ {
2136	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
2137	return it2_.index1 ();
2138	}
2139	}
2140	BOOST_UBLAS_INLINE
2141	size_type index2 () const {
2142	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2143	if (current_ == 0) {
2144	BOOST_UBLAS_CHECK (it1_ != it1_end_, internal_logic ());
2145	return it1_.index1 ();
2146	} else /* if (current_ == 1) */ {
2147	BOOST_UBLAS_CHECK (it2_ != it2_end_, internal_logic ());
2148	return it2_.index2 ();
2149	}
2150	}
2151
2152	// Assignment
2153	BOOST_UBLAS_INLINE
2154	const_iterator2 &operator = (const const_iterator2 &it) {
2155	container_const_reference<self_type>::assign (&it ());
2156	begin_ = it.begin_;
2157	end_ = it.end_;
2158	current_ = it.current_;
2159	it1_begin_ = it.it1_begin_;
2160	it1_end_ = it.it1_end_;
2161	it1_ = it.it1_;
2162	it2_begin_ = it.it2_begin_;
2163	it2_end_ = it.it2_end_;
2164	it2_ = it.it2_;
2165	return *this;
2166	}
2167
2168	// Comparison
2169	BOOST_UBLAS_INLINE
2170	bool operator == (const const_iterator2 &it) const {
2171	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2172	BOOST_UBLAS_CHECK (current_ == 0 \|\| current_ == 1, internal_logic ());
2173	BOOST_UBLAS_CHECK (it.current_ == 0 \|\| it.current_ == 1, internal_logic ());
2174	BOOST_UBLAS_CHECK (/* begin_ == it.begin_ && */ end_ == it.end_, internal_logic ());
2175	return (current_ == 0 && it.current_ == 0 && it1_ == it.it1_) \|\|
2176	(current_ == 1 && it.current_ == 1 && it2_ == it.it2_);
2177	}
2178	BOOST_UBLAS_INLINE
2179	bool operator < (const const_iterator2 &it) const {
2180	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2181	return it - *this > 0;
2182	}
2183
2184	private:
2185	int begin_;
2186	int end_;
2187	int current_;
2188	const_subiterator1_type it1_begin_;
2189	const_subiterator1_type it1_end_;
2190	const_subiterator1_type it1_;
2191	const_subiterator2_type it2_begin_;
2192	const_subiterator2_type it2_end_;
2193	const_subiterator2_type it2_;
2194	};
2195	#endif
2196
2197	BOOST_UBLAS_INLINE
2198	const_iterator2 begin2 () const {
2199	return find2 (0, 0, 0);
2200	}
2201	BOOST_UBLAS_INLINE
2202	const_iterator2 end2 () const {
2203	return find2 (0, 0, size2 ());
2204	}
2205
2206	#ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2207	class iterator2:
2208	public container_reference<hermitian_adaptor>,
2209	public random_access_iterator_base<typename iterator_restrict_traits<
2210	typename subiterator2_type::iterator_category, packed_random_access_iterator_tag>::iterator_category,
2211	iterator2, value_type> {
2212	public:
2213	typedef typename subiterator2_type::value_type value_type;
2214	typedef typename subiterator2_type::difference_type difference_type;
2215	typedef typename subiterator2_type::reference reference;
2216	typedef typename subiterator2_type::pointer pointer;
2217
2218	typedef iterator1 dual_iterator_type;
2219	typedef reverse_iterator1 dual_reverse_iterator_type;
2220
2221	// Construction and destruction
2222	BOOST_UBLAS_INLINE
2223	iterator2 ():
2224	container_reference<self_type> (), it2_ () {}
2225	BOOST_UBLAS_INLINE
2226	iterator2 (self_type &m, const subiterator2_type &it2):
2227	container_reference<self_type> (m), it2_ (it2) {}
2228
2229	// Arithmetic
2230	BOOST_UBLAS_INLINE
2231	iterator2 &operator ++ () {
2232	++ it2_;
2233	return *this;
2234	}
2235	BOOST_UBLAS_INLINE
2236	iterator2 &operator -- () {
2237	-- it2_;
2238	return *this;
2239	}
2240	BOOST_UBLAS_INLINE
2241	iterator2 &operator += (difference_type n) {
2242	it2_ += n;
2243	return *this;
2244	}
2245	BOOST_UBLAS_INLINE
2246	iterator2 &operator -= (difference_type n) {
2247	it2_ -= n;
2248	return *this;
2249	}
2250	BOOST_UBLAS_INLINE
2251	difference_type operator - (const iterator2 &it) const {
2252	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2253	return it2_ - it.it2_;
2254	}
2255
2256	// Dereference
2257	BOOST_UBLAS_INLINE
2258	reference operator * () const {
2259	return *it2_;
2260	}
2261
2262	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2263	BOOST_UBLAS_INLINE
2264	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2265	typename self_type::
2266	#endif
2267	iterator1 begin () const {
2268	return (*this) ().find1 (1, 0, index2 ());
2269	}
2270	BOOST_UBLAS_INLINE
2271	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2272	typename self_type::
2273	#endif
2274	iterator1 end () const {
2275	return (this) ().find1 (1, (this) ().size1 (), index2 ());
2276	}
2277	BOOST_UBLAS_INLINE
2278	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2279	typename self_type::
2280	#endif
2281	reverse_iterator1 rbegin () const {
2282	return reverse_iterator1 (end ());
2283	}
2284	BOOST_UBLAS_INLINE
2285	#ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2286	typename self_type::
2287	#endif
2288	reverse_iterator1 rend () const {
2289	return reverse_iterator1 (begin ());
2290	}
2291	#endif
2292
2293	// Indices
2294	BOOST_UBLAS_INLINE
2295	size_type index1 () const {
2296	return it2_.index1 ();
2297	}
2298	BOOST_UBLAS_INLINE
2299	size_type index2 () const {
2300	return it2_.index2 ();
2301	}
2302
2303	// Assignment
2304	BOOST_UBLAS_INLINE
2305	iterator2 &operator = (const iterator2 &it) {
2306	container_reference<self_type>::assign (&it ());
2307	it2_ = it.it2_;
2308	return *this;
2309	}
2310
2311	// Comparison
2312	BOOST_UBLAS_INLINE
2313	bool operator == (const iterator2 &it) const {
2314	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2315	return it2_ == it.it2_;
2316	}
2317	BOOST_UBLAS_INLINE
2318	bool operator < (const iterator2 &it) const {
2319	BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2320	return it2_ < it.it2_;
2321	}
2322
2323	private:
2324	subiterator2_type it2_;
2325
2326	friend class const_iterator2;
2327	};
2328	#endif
2329
2330	BOOST_UBLAS_INLINE
2331	iterator2 begin2 () {
2332	return find2 (0, 0, 0);
2333	}
2334	BOOST_UBLAS_INLINE
2335	iterator2 end2 () {
2336	return find2 (0, 0, size2 ());
2337	}
2338
2339	// Reverse iterators
2340
2341	BOOST_UBLAS_INLINE
2342	const_reverse_iterator1 rbegin1 () const {
2343	return const_reverse_iterator1 (end1 ());
2344	}
2345	BOOST_UBLAS_INLINE
2346	const_reverse_iterator1 rend1 () const {
2347	return const_reverse_iterator1 (begin1 ());
2348	}
2349
2350	BOOST_UBLAS_INLINE
2351	reverse_iterator1 rbegin1 () {
2352	return reverse_iterator1 (end1 ());
2353	}
2354	BOOST_UBLAS_INLINE
2355	reverse_iterator1 rend1 () {
2356	return reverse_iterator1 (begin1 ());
2357	}
2358
2359	BOOST_UBLAS_INLINE
2360	const_reverse_iterator2 rbegin2 () const {
2361	return const_reverse_iterator2 (end2 ());
2362	}
2363	BOOST_UBLAS_INLINE
2364	const_reverse_iterator2 rend2 () const {
2365	return const_reverse_iterator2 (begin2 ());
2366	}
2367
2368	BOOST_UBLAS_INLINE
2369	reverse_iterator2 rbegin2 () {
2370	return reverse_iterator2 (end2 ());
2371	}
2372	BOOST_UBLAS_INLINE
2373	reverse_iterator2 rend2 () {
2374	return reverse_iterator2 (begin2 ());
2375	}
2376
2377	private:
2378	matrix_closure_type data_;
2379	static value_type conj_;
2380	};
2381
2382	template<class M, class TRI>
2383	typename hermitian_adaptor<M, TRI>::value_type hermitian_adaptor<M, TRI>::conj_;
2384
2385	// Specialization for temporary_traits
2386	template <class M, class TRI>
2387	struct vector_temporary_traits< hermitian_adaptor<M, TRI> >
2388	: vector_temporary_traits< M > {} ;
2389
2390	template <class M, class TRI>
2391	struct matrix_temporary_traits< hermitian_adaptor<M, TRI> >
2392	: matrix_temporary_traits< M > {} ;
2393
2394	}}}
2395
2396	#endif

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: NonGTP/Boost/boost/numeric/ublas/hermitian.hpp @ 857

Download in other formats: