Blaze 3.9
SymmetricMatrix.h
Go to the documentation of this file.
1//=================================================================================================
33//=================================================================================================
34
35#ifndef _BLAZE_MATH_SYMMETRICMATRIX_H_
36#define _BLAZE_MATH_SYMMETRICMATRIX_H_
37
38
39//*************************************************************************************************
40// Includes
41//*************************************************************************************************
42
43#include <cmath>
44#include <vector>
45#include <blaze/math/Aliases.h>
58#include <blaze/util/Assert.h>
60#include <blaze/util/Random.h>
61#include <blaze/util/Types.h>
62
63
64namespace blaze {
65
66//=================================================================================================
67//
68// RAND SPECIALIZATION
69//
70//=================================================================================================
71
72//*************************************************************************************************
79template< typename MT // Type of the adapted matrix
80 , bool SO // Storage order of the adapted matrix
81 , bool DF // Density flag
82 , bool SF > // Scalar flag
83class Rand< SymmetricMatrix<MT,SO,DF,SF> >
84{
85 public:
86 //**********************************************************************************************
91 inline const SymmetricMatrix<MT,SO,DF,SF> generate() const
92 {
94
95 SymmetricMatrix<MT,SO,DF,SF> matrix;
96 randomize( matrix );
97 return matrix;
98 }
99 //**********************************************************************************************
100
101 //**********************************************************************************************
107 inline const SymmetricMatrix<MT,SO,DF,SF> generate( size_t n ) const
108 {
110
111 SymmetricMatrix<MT,SO,DF,SF> matrix( n );
112 randomize( matrix );
113 return matrix;
114 }
115 //**********************************************************************************************
116
117 //**********************************************************************************************
125 inline const SymmetricMatrix<MT,SO,DF,SF> generate( size_t n, size_t nonzeros ) const
126 {
129
130 if( nonzeros > n*n ) {
131 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
132 }
133
134 SymmetricMatrix<MT,SO,DF,SF> matrix( n );
135 randomize( matrix, nonzeros );
136
137 return matrix;
138 }
139 //**********************************************************************************************
140
141 //**********************************************************************************************
148 template< typename Arg > // Min/max argument type
149 inline const SymmetricMatrix<MT,SO,DF,SF> generate( const Arg& min, const Arg& max ) const
150 {
152
153 SymmetricMatrix<MT,SO,DF,SF> matrix;
154 randomize( matrix, min, max );
155 return matrix;
156 }
157 //**********************************************************************************************
158
159 //**********************************************************************************************
167 template< typename Arg > // Min/max argument type
168 inline const SymmetricMatrix<MT,SO,DF,SF> generate( size_t n, const Arg& min, const Arg& max ) const
169 {
171
172 SymmetricMatrix<MT,SO,DF,SF> matrix( n );
173 randomize( matrix, min, max );
174 return matrix;
175 }
176 //**********************************************************************************************
177
178 //**********************************************************************************************
188 template< typename Arg > // Min/max argument type
189 inline const SymmetricMatrix<MT,SO,DF,SF>
190 generate( size_t n, size_t nonzeros, const Arg& min, const Arg& max ) const
191 {
194
195 if( nonzeros > n*n ) {
196 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
197 }
198
199 SymmetricMatrix<MT,SO,DF,SF> matrix( n );
200 randomize( matrix, nonzeros, min, max );
201
202 return matrix;
203 }
204 //**********************************************************************************************
205
206 //**********************************************************************************************
212 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix ) const
213 {
214 randomize( matrix, typename IsDenseMatrix<MT>::Type() );
215 }
216 //**********************************************************************************************
217
218 //**********************************************************************************************
226 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix, size_t nonzeros ) const
227 {
229
230 using ET = ElementType_t<MT>;
231
232 const size_t n( matrix.rows() );
233
234 if( nonzeros > n*n ) {
235 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
236 }
237
238 if( n == 0UL ) return;
239
240 matrix.reset();
241 matrix.reserve( nonzeros+1UL );
242
243 while( matrix.nonZeros() < nonzeros ) {
244 matrix( rand<size_t>( 0UL, n-1UL ), rand<size_t>( 0UL, n-1UL ) ) = rand<ET>();
245 }
246 }
247 //**********************************************************************************************
248
249 //**********************************************************************************************
257 template< typename Arg > // Min/max argument type
258 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix,
259 const Arg& min, const Arg& max ) const
260 {
261 randomize( matrix, min, max, typename IsDenseMatrix<MT>::Type() );
262 }
263 //**********************************************************************************************
264
265 //**********************************************************************************************
275 template< typename Arg > // Min/max argument type
276 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix,
277 size_t nonzeros, const Arg& min, const Arg& max ) const
278 {
280
281 using ET = ElementType_t<MT>;
282
283 const size_t n( matrix.rows() );
284
285 if( nonzeros > n*n ) {
286 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
287 }
288
289 if( n == 0UL ) return;
290
291 std::vector<size_t> dist( n );
292 std::vector<bool> structure( n*n );
293 size_t nz( 0UL );
294
295 while( nz < nonzeros )
296 {
297 const size_t row = rand<size_t>( 0UL, n-1UL );
298 const size_t col = rand<size_t>( 0UL, n-1UL );
299
300 if( structure[row*n+col] ) continue;
301
302 ++dist[row];
303 structure[row*n+col] = true;
304 ++nz;
305
306 if( row != col ) {
307 ++dist[col];
308 structure[col*n+row] = true;
309 ++nz;
310 }
311 }
312
313 matrix.reset();
314 matrix.reserve( nz );
315
316 for( size_t i=0UL; i<n; ++i ) {
317 matrix.reserve( i, dist[i] );
318 }
319
320 for( size_t i=0UL; i<n; ++i ) {
321 for( size_t j=i; j<n; ++j ) {
322 if( structure[i*n+j] ) {
323 matrix.append( i, j, rand<ET>( min, max ) );
324 }
325 }
326 }
327 }
328 //**********************************************************************************************
329
330 private:
331 //*************************************************************************************************
337 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix, TrueType ) const
338 {
340
341 using ET = ElementType_t<MT>;
342
343 const size_t n( matrix.rows() );
344
345 for( size_t i=0UL; i<n; ++i ) {
346 for( size_t j=0UL; j<=i; ++j ) {
347 matrix(i,j) = rand<ET>();
348 }
349 }
350 }
351 //*************************************************************************************************
352
353 //*************************************************************************************************
359 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix, FalseType ) const
360 {
362
363 const size_t n( matrix.rows() );
364
365 if( n == 0UL ) return;
366
367 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.5*n*n ) ) );
368
369 randomize( matrix, nonzeros );
370 }
371 //*************************************************************************************************
372
373 //*************************************************************************************************
381 template< typename Arg > // Min/max argument type
382 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix,
383 const Arg& min, const Arg& max, TrueType ) const
384 {
386
387 using ET = ElementType_t<MT>;
388
389 const size_t n( matrix.rows() );
390
391 for( size_t i=0UL; i<n; ++i ) {
392 for( size_t j=0UL; j<=i; ++j ) {
393 matrix(i,j) = rand<ET>( min, max );
394 }
395 }
396 }
397 //*************************************************************************************************
398
399 //*************************************************************************************************
407 template< typename Arg > // Min/max argument type
408 inline void randomize( SymmetricMatrix<MT,SO,DF,SF>& matrix,
409 const Arg& min, const Arg& max, FalseType ) const
410 {
412
413 const size_t n( matrix.rows() );
414
415 if( n == 0UL ) return;
416
417 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.5*n*n ) ) );
418
419 randomize( matrix, nonzeros, min, max );
420 }
421 //*************************************************************************************************
422};
424//*************************************************************************************************
425
426
427
428
429//=================================================================================================
430//
431// MAKE FUNCTIONS
432//
433//=================================================================================================
434
435//*************************************************************************************************
442template< typename MT // Type of the adapted matrix
443 , bool SO // Storage order of the adapted matrix
444 , bool SF > // Scalar flag
445void makeSymmetric( SymmetricMatrix<MT,SO,true,SF>& matrix )
446{
447 using blaze::randomize;
448
449 randomize( matrix );
450}
452//*************************************************************************************************
453
454
455//*************************************************************************************************
464template< typename MT // Type of the adapted matrix
465 , bool SO // Storage order of the adapted matrix
466 , bool SF // Scalar flag
467 , typename Arg > // Min/max argument type
468void makeSymmetric( SymmetricMatrix<MT,SO,true,SF>& matrix, const Arg& min, const Arg& max )
469{
470 using blaze::randomize;
471
472 randomize( matrix, min, max );
473}
475//*************************************************************************************************
476
477
478//*************************************************************************************************
485template< typename MT // Type of the adapted matrix
486 , bool SO // Storage order of the adapted matrix
487 , bool SF > // Scalar flag
488void makeHermitian( SymmetricMatrix<MT,SO,true,SF>& matrix )
489{
490 using BT = UnderlyingBuiltin_t< ElementType_t<MT> >;
491
493
494 const size_t n( matrix.rows() );
495
496 for( size_t i=0UL; i<n; ++i ) {
497 for( size_t j=0UL; j<=i; ++j ) {
498 matrix(i,j) = rand<BT>();
499 }
500 }
501
502 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
503}
505//*************************************************************************************************
506
507
508//*************************************************************************************************
517template< typename MT // Type of the adapted matrix
518 , bool SO // Storage order of the adapted matrix
519 , bool SF // Scalar flag
520 , typename Arg > // Min/max argument type
521void makeHermitian( SymmetricMatrix<MT,SO,true,SF>& matrix, const Arg& min, const Arg& max )
522{
523 using BT = UnderlyingBuiltin_t< ElementType_t<MT> >;
524
526
527 const size_t n( matrix.rows() );
528
529 for( size_t i=0UL; i<n; ++i ) {
530 for( size_t j=0UL; j<=i; ++j ) {
531 matrix(i,j) = rand<BT>( real( min ), real( max ) );
532 }
533 }
534
535 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
536}
538//*************************************************************************************************
539
540
541//*************************************************************************************************
548template< typename MT // Type of the adapted matrix
549 , bool SO // Storage order of the adapted matrix
550 , bool SF > // Scalar flag
551void makePositiveDefinite( SymmetricMatrix<MT,SO,true,SF>& matrix )
552{
553 using BT = UnderlyingBuiltin_t< ElementType_t<MT> >;
554
556
557 const size_t n( matrix.rows() );
558
559 makeHermitian( matrix );
560 matrix *= matrix;
561
562 for( size_t i=0UL; i<n; ++i ) {
563 matrix(i,i) += BT(n);
564 }
565
566 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
567}
569//*************************************************************************************************
570
571} // namespace blaze
572
573#endif
Header file for auxiliary alias declarations.
Header file for run time assertion macros.
Header file for all basic DenseMatrix functionality.
Header file for the IntegralConstant class template.
Header file for the IsDenseMatrix type trait.
Constraint on the data type.
Constraint on the data type.
Header file for all basic SparseMatrix functionality.
Header file for the UnderlyingBuiltin type trait.
Header file for the implementation of a diagonal matrix adaptor.
Header file for the implementation of a symmetric matrix adaptor.
Constraint on the data type.
Constraint on the data type.
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:1339
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:1375
decltype(auto) real(const DenseMatrix< MT, SO > &dm)
Returns a matrix containing the real part of each single element of dm.
Definition: DMatMapExpr.h:1529
decltype(auto) ceil(const DenseMatrix< MT, SO > &dm)
Applies the ceil() function to each single element of the dense matrix dm.
Definition: DMatMapExpr.h:1380
bool isHermitian(const DenseMatrix< MT, SO > &dm)
Checks if the given dense matrix is Hermitian.
Definition: DenseMatrix.h:1534
decltype(auto) generate(size_t m, size_t n, OP op)
Generates a new dense matrix filled via the given custom binary operation.
Definition: DMatGenExpr.h:675
#define BLAZE_CONSTRAINT_MUST_BE_DENSE_MATRIX_TYPE(T)
Constraint on the data type.
Definition: DenseMatrix.h:61
#define BLAZE_CONSTRAINT_MUST_NOT_BE_RESIZABLE_TYPE(T)
Constraint on the data type.
Definition: Resizable.h:81
#define BLAZE_CONSTRAINT_MUST_BE_SPARSE_MATRIX_TYPE(T)
Constraint on the data type.
Definition: SparseMatrix.h:61
#define BLAZE_CONSTRAINT_MUST_BE_RESIZABLE_TYPE(T)
Constraint on the data type.
Definition: Resizable.h:61
#define BLAZE_CONSTRAINT_MUST_BE_SCALAR_TYPE(T)
Constraint on the data type.
Definition: Scalar.h:61
void randomize(T &&value)
Randomization of a given variable.
Definition: Random.h:626
decltype(auto) row(Matrix< MT, SO > &, RRAs...)
Creating a view on a specific row of the given matrix.
Definition: Row.h:137
#define BLAZE_INTERNAL_ASSERT(expr, msg)
Run time assertion macro for internal checks.
Definition: Assert.h:101
BoolConstant< true > TrueType
Type traits base class.
Definition: IntegralConstant.h:132
BoolConstant< false > FalseType
Type/value traits base class.
Definition: IntegralConstant.h:121
#define BLAZE_THROW_INVALID_ARGUMENT(MESSAGE)
Macro for the emission of a std::invalid_argument exception.
Definition: Exception.h:235
Header file for the exception macros of the math module.
Header file for the real shim.
Implementation of a random number generator.
Header file for basic type definitions.