CustomMatrix.h
Go to the documentation of this file.
1 //=================================================================================================
33 //=================================================================================================
34 
35 #ifndef _BLAZE_MATH_CUSTOMMATRIX_H_
36 #define _BLAZE_MATH_CUSTOMMATRIX_H_
37 
38 
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42 
45 #include <blaze/math/DenseMatrix.h>
47 #include <blaze/math/Exception.h>
50 #include <blaze/math/shims/Real.h>
52 #include <blaze/math/ZeroMatrix.h>
53 #include <blaze/util/Assert.h>
55 #include <blaze/util/Random.h>
56 
57 
58 namespace blaze {
59 
60 //=================================================================================================
61 //
62 // RAND SPECIALIZATION
63 //
64 //=================================================================================================
65 
66 //*************************************************************************************************
73 template< typename Type // Data type of the matrix
74  , bool AF // Alignment flag
75  , bool PF // Padding flag
76  , bool SO // Storage order
77  , typename RT > // Result type
78 class Rand< CustomMatrix<Type,AF,PF,SO,RT> >
79 {
80  public:
81  //**Randomize functions*************************************************************************
84  inline void randomize( CustomMatrix<Type,AF,PF,SO,RT>& matrix ) const;
85 
86  template< typename Arg >
87  inline void randomize( CustomMatrix<Type,AF,PF,SO,RT>& matrix, const Arg& min, const Arg& max ) const;
89  //**********************************************************************************************
90 };
92 //*************************************************************************************************
93 
94 
95 //*************************************************************************************************
102 template< typename Type // Data type of the matrix
103  , bool AF // Alignment flag
104  , bool PF // Padding flag
105  , bool SO // Storage order
106  , typename RT > // Result type
107 inline void Rand< CustomMatrix<Type,AF,PF,SO,RT> >::randomize( CustomMatrix<Type,AF,PF,SO,RT>& matrix ) const
108 {
109  using blaze::randomize;
110 
111  const size_t m( matrix.rows() );
112  const size_t n( matrix.columns() );
113 
114  for( size_t i=0UL; i<m; ++i ) {
115  for( size_t j=0UL; j<n; ++j ) {
116  randomize( matrix(i,j) );
117  }
118  }
119 }
121 //*************************************************************************************************
122 
123 
124 //*************************************************************************************************
133 template< typename Type // Data type of the matrix
134  , bool AF // Alignment flag
135  , bool PF // Padding flag
136  , bool SO // Storage order
137  , typename RT > // Result type
138 template< typename Arg > // Min/max argument type
139 inline void Rand< CustomMatrix<Type,AF,PF,SO,RT> >::randomize( CustomMatrix<Type,AF,PF,SO,RT>& matrix,
140  const Arg& min, const Arg& max ) const
141 {
142  using blaze::randomize;
143 
144  const size_t m( matrix.rows() );
145  const size_t n( matrix.columns() );
146 
147  for( size_t i=0UL; i<m; ++i ) {
148  for( size_t j=0UL; j<n; ++j ) {
149  randomize( matrix(i,j), min, max );
150  }
151  }
152 }
154 //*************************************************************************************************
155 
156 
157 
158 
159 //=================================================================================================
160 //
161 // MAKE FUNCTIONS
162 //
163 //=================================================================================================
164 
165 //*************************************************************************************************
173 template< typename Type // Data type of the matrix
174  , bool AF // Alignment flag
175  , bool PF // Padding flag
176  , bool SO // Storage order
177  , typename RT > // Result type
178 void makeSymmetric( CustomMatrix<Type,AF,PF,SO,RT>& matrix )
179 {
180  using blaze::randomize;
181 
182  if( !isSquare( ~matrix ) ) {
183  BLAZE_THROW_INVALID_ARGUMENT( "Invalid non-square matrix provided" );
184  }
185 
186  const size_t n( matrix.rows() );
187 
188  for( size_t i=0UL; i<n; ++i ) {
189  for( size_t j=0UL; j<i; ++j ) {
190  randomize( matrix(i,j) );
191  matrix(j,i) = matrix(i,j);
192  }
193  randomize( matrix(i,i) );
194  }
195 
196  BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
197 }
199 //*************************************************************************************************
200 
201 
202 //*************************************************************************************************
212 template< typename Type // Data type of the matrix
213  , bool AF // Alignment flag
214  , bool PF // Padding flag
215  , bool SO // Storage order
216  , typename RT // Result type
217  , typename Arg > // Min/max argument type
218 void makeSymmetric( CustomMatrix<Type,AF,PF,SO,RT>& matrix, const Arg& min, const Arg& max )
219 {
220  using blaze::randomize;
221 
222  if( !isSquare( ~matrix ) ) {
223  BLAZE_THROW_INVALID_ARGUMENT( "Invalid non-square matrix provided" );
224  }
225 
226  const size_t n( matrix.rows() );
227 
228  for( size_t i=0UL; i<n; ++i ) {
229  for( size_t j=0UL; j<i; ++j ) {
230  randomize( matrix(i,j), min, max );
231  matrix(j,i) = matrix(i,j);
232  }
233  randomize( matrix(i,i), min, max );
234  }
235 
236  BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
237 }
239 //*************************************************************************************************
240 
241 
242 //*************************************************************************************************
250 template< typename Type // Data type of the matrix
251  , bool AF // Alignment flag
252  , bool PF // Padding flag
253  , bool SO // Storage order
254  , typename RT > // Result type
255 void makeHermitian( CustomMatrix<Type,AF,PF,SO,RT>& matrix )
256 {
257  using blaze::randomize;
258 
260 
261  using BT = UnderlyingBuiltin_t<Type>;
262 
263  if( !isSquare( ~matrix ) ) {
264  BLAZE_THROW_INVALID_ARGUMENT( "Invalid non-square matrix provided" );
265  }
266 
267  const size_t n( matrix.rows() );
268 
269  for( size_t i=0UL; i<n; ++i ) {
270  for( size_t j=0UL; j<i; ++j ) {
271  randomize( matrix(i,j) );
272  matrix(j,i) = conj( matrix(i,j) );
273  }
274  matrix(i,i) = rand<BT>();
275  }
276 
277  BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
278 }
280 //*************************************************************************************************
281 
282 
283 //*************************************************************************************************
293 template< typename Type // Data type of the matrix
294  , bool AF // Alignment flag
295  , bool PF // Padding flag
296  , bool SO // Storage order
297  , typename RT // Result type
298  , typename Arg > // Min/max argument type
299 void makeHermitian( CustomMatrix<Type,AF,PF,SO,RT>& matrix, const Arg& min, const Arg& max )
300 {
301  using blaze::randomize;
302 
304 
305  using BT = UnderlyingBuiltin_t<Type>;
306 
307  if( !isSquare( ~matrix ) ) {
308  BLAZE_THROW_INVALID_ARGUMENT( "Invalid non-square matrix provided" );
309  }
310 
311  const size_t n( matrix.rows() );
312 
313  for( size_t i=0UL; i<n; ++i ) {
314  for( size_t j=0UL; j<i; ++j ) {
315  randomize( matrix(i,j), min, max );
316  matrix(j,i) = conj( matrix(i,j) );
317  }
318  matrix(i,i) = rand<BT>( real( min ), real( max ) );
319  }
320 
321  BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
322 }
324 //*************************************************************************************************
325 
326 
327 //*************************************************************************************************
335 template< typename Type // Data type of the matrix
336  , bool AF // Alignment flag
337  , bool PF // Padding flag
338  , bool SO // Storage order
339  , typename RT > // Result type
340 void makePositiveDefinite( CustomMatrix<Type,AF,PF,SO,RT>& matrix )
341 {
342  using blaze::randomize;
343 
345 
346  if( !isSquare( ~matrix ) ) {
347  BLAZE_THROW_INVALID_ARGUMENT( "Invalid non-square matrix provided" );
348  }
349 
350  const size_t n( matrix.rows() );
351 
352  randomize( matrix );
353  matrix *= ctrans( matrix );
354 
355  for( size_t i=0UL; i<n; ++i ) {
356  matrix(i,i) += Type(n);
357  }
358 
359  BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-symmetric matrix detected" );
360 }
362 //*************************************************************************************************
363 
364 } // namespace blaze
365 
366 #endif
#define BLAZE_THROW_INVALID_ARGUMENT(MESSAGE)
Macro for the emission of a std::invalid_argument exception.This macro encapsulates the default way o...
Definition: Exception.h:235
Constraint on the data type.
decltype(auto) real(const DenseMatrix< MT, SO > &dm)
Returns a matrix containing the real part of each single element of dm.
Definition: DMatMapExpr.h:1392
void randomize(T &&value)
Randomization of a given variable.
Definition: Random.h:929
Implementation of a random number generator.
Header file for the implementation of a customizable matrix.
void randomize(T &value) const
Randomization of the given variable with a new value in the range .
Definition: Random.h:292
Namespace of the Blaze C++ math library.
Definition: Blaze.h:58
decltype(auto) ctrans(const DenseMatrix< MT, SO > &dm)
Returns the conjugate transpose matrix of dm.
Definition: DMatMapExpr.h:1364
Header file for the UnderlyingBuiltin type trait.
decltype(auto) min(const DenseMatrix< MT1, SO1 > &lhs, const DenseMatrix< MT2, SO2 > &rhs)
Computes the componentwise minimum of the dense matrices lhs and rhs.
Definition: DMatDMatMapExpr.h:1147
Header file for the exception macros of the math module.
Header file for the implementation of a dynamic MxN matrix.
decltype(auto) max(const DenseMatrix< MT1, SO1 > &lhs, const DenseMatrix< MT2, SO2 > &rhs)
Computes the componentwise maximum of the dense matrices lhs and rhs.
Definition: DMatDMatMapExpr.h:1179
Header file for the complete DynamicVector implementation.
Header file for the conjugate shim.
Header file for run time assertion macros.
#define BLAZE_CONSTRAINT_MUST_BE_NUMERIC_TYPE(T)
Constraint on the data type.In case the given data type T is not a numeric (integral or floating poin...
Definition: Numeric.h:61
bool isHermitian(const DenseMatrix< MT, SO > &dm)
Checks if the given dense matrix is Hermitian.
Definition: DenseMatrix.h:617
bool isSymmetric(const DenseMatrix< MT, SO > &dm)
Checks if the given dense matrix is symmetric.
Definition: DenseMatrix.h:539
Header file for all basic DenseMatrix functionality.
Header file for the complete IdentityMatrix implementation.
Header file for the real shim.
decltype(auto) conj(const DenseMatrix< MT, SO > &dm)
Returns a matrix containing the complex conjugate of each single element of dm.
Definition: DMatMapExpr.h:1326
bool isSquare(const Matrix< MT, SO > &matrix) noexcept
Checks if the given matrix is a square matrix.
Definition: Matrix.h:951
#define BLAZE_INTERNAL_ASSERT(expr, msg)
Run time assertion macro for internal checks.In case of an invalid run time expression, the program execution is terminated. The BLAZE_INTERNAL_ASSERT macro can be disabled by setting the BLAZE_USER_ASSERTION flag to zero or by defining NDEBUG during the compilation.
Definition: Assert.h:101
Header file for the complete ZeroMatrix implementation.