1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
|
/*
*
* Copyright (c) Karl Meerbergen & Kresimir Fresl 2003
*
* Distributed under the Boost Software License, Version 1.0.
* (See accompanying file LICENSE_1_0.txt or copy at
* http://www.boost.org/LICENSE_1_0.txt)
*
* KF acknowledges the support of the Faculty of Civil Engineering,
* University of Zagreb, Croatia.
*
*/
#ifndef BOOST_NUMERIC_BINDINGS_LAPACK_WORKSPACE_HPP
#define BOOST_NUMERIC_BINDINGS_LAPACK_WORKSPACE_HPP
#include <boost/mpl/bool.hpp>
#include <boost/numeric/bindings/value_type.hpp>
#include <memory>
namespace boost {
namespace numeric {
namespace bindings {
namespace lapack {
/*
* Organization of workspace in Lapack.
* We allow one of the following arguments in a number of Lapack functions
* - minimal_workspace() : the function allocates the minimum workspace required for the function
* - optimal_workspace() : the function allocates the amount of workspace that allows optimal
* execution.
* - workspace( work ) : the function uses the workspace array in work.
* - workspace( rwork, work ) : the function uses a real array rwork and a compolex array work as
* workspace. (There are Lapack functions for complex matrices
* that require two workarrays)
* */
// Four classes are introduced to distinguish between the different type of memory allocations
struct minimal_workspace {} ;
struct optimal_workspace {} ;
namespace detail {
template <typename W>
class workspace1 {
public:
workspace1(W& w)
: w_( w )
{}
public:
typedef typename bindings::value_type< W>::type value_type ;
W& select( value_type const& ) { return w_ ; }
private:
W& w_ ;
}; // struct workspace1
template <typename W, typename WRI>
class workspace2 {
public:
workspace2(W& w, WRI& wri)
: w_(w)
, wri_(wri)
{}
public:
typedef typename bindings::value_type< W>::type w_value_type ;
W& select( w_value_type const& ) { return w_ ; }
typedef typename bindings::value_type< WRI>::type wri_value_type ;
WRI& select( wri_value_type const& ) { return wri_ ; }
private:
W& w_ ;
WRI& wri_ ;
}; // struct workspace2
template <typename W, typename WR, typename WI>
class workspace3 {
public:
workspace3(W& w, WR& wr, WI& wi)
: w_(w)
, wr_(wr)
, wi_(wi)
{}
public:
typedef typename bindings::value_type< W>::type w_value_type ;
W& select( w_value_type const& ) { return w_ ; }
typedef typename bindings::value_type< WR>::type wr_value_type ;
WR& select( wr_value_type const& ) { return wr_ ; }
typedef typename bindings::value_type< WI>::type wi_value_type ;
WI& select( wi_value_type const& ) { return wi_ ; }
private:
W& w_ ;
WR& wr_ ;
WI& wi_ ;
}; // struct workspace3
template <typename W, typename WR, typename WI, typename WB>
class workspace4 {
public:
workspace4(W& w, WR& wr, WI& wi, WB& wb)
: w_(w)
, wr_(wr)
, wi_(wi)
, wb_(wb)
{}
public:
typedef typename bindings::value_type< W>::type w_value_type ;
W& select( w_value_type const& ) { return w_ ; }
typedef typename bindings::value_type< WR>::type wr_value_type ;
WR& select( wr_value_type const& ) { return wr_ ; }
typedef typename bindings::value_type< WI>::type wi_value_type ;
WI& select( wi_value_type const& ) { return wi_ ; }
typedef typename bindings::value_type< WB>::type wb_value_type ;
WB& select( wb_value_type const& ) { return wb_ ; }
private:
W& w_ ;
WR& wr_ ;
WI& wi_ ;
WB& wb_ ;
}; // struct workspace4
}
template <typename W>
detail::workspace1<W> workspace(W& w) {
return detail::workspace1<W>(w) ;
} // workspace()
//
// Two situations:
// Real valued: workspace( real array, integer array )
// Complex valued: workspace( complex array, real array )
//
template <typename W, typename WRI>
detail::workspace2<W,WRI> workspace(W& w, WRI& wri) {
return detail::workspace2<W,WRI>(w, wri) ;
} // workspace()
//
// Complex valued: workspace( complex array, real array, integer array )
//
template <typename W, typename WR, typename WI>
detail::workspace3<W,WR,WI> workspace(W& w, WR& wr, WI& wi) {
return detail::workspace3<W,WR,WI>(w, wr, wi) ;
} // workspace()
//
// Complex valued: workspace( complex array, real array, integer array, bool array )
//
template <typename W, typename WR, typename WI, typename WB>
detail::workspace4<W,WR,WI,WB> workspace(W& w, WR& wr, WI& wi, WB& wb) {
return detail::workspace4<W,WR,WI,WB>(w, wr, wi, wb) ;
} // workspace()
namespace detail {
template< typename T >
struct is_workspace: mpl::false_ {};
template<>
struct is_workspace< minimal_workspace >: mpl::true_ {};
template<>
struct is_workspace< optimal_workspace >: mpl::true_ {};
template< typename T >
struct is_workspace< workspace1<T> >: mpl::true_ {};
template< typename T1, typename T2 >
struct is_workspace< workspace2<T1,T2> >: mpl::true_ {};
template< typename T1, typename T2, typename T3 >
struct is_workspace< workspace3<T1,T2,T3> >: mpl::true_ {};
template< typename T1, typename T2, typename T3, typename T4 >
struct is_workspace< workspace4<T1,T2,T3,T4> >: mpl::true_ {};
}
} // namespace lapack
} // namespace bindings
} // namespace numeric
} // namespace boost
#endif
|
