Blaze 3.9
UpperMatrix.h
Go to the documentation of this file.
1//=================================================================================================
33//=================================================================================================
34
35#ifndef _BLAZE_MATH_UPPERMATRIX_H_
36#define _BLAZE_MATH_UPPERMATRIX_H_
37
38
39//*************************************************************************************************
40// Includes
41//*************************************************************************************************
42
43#include <cmath>
44#include <vector>
45#include <blaze/math/Aliases.h>
58#include <blaze/util/Indices.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
82class Rand< UpperMatrix<MT,SO,DF> >
83{
84 public:
85 //**********************************************************************************************
90 inline const UpperMatrix<MT,SO,DF> generate() const
91 {
93
94 UpperMatrix<MT,SO,DF> matrix;
95 randomize( matrix );
96 return matrix;
97 }
98 //**********************************************************************************************
99
100 //**********************************************************************************************
106 inline const UpperMatrix<MT,SO,DF> generate( size_t n ) const
107 {
109
110 UpperMatrix<MT,SO,DF> matrix( n );
111 randomize( matrix );
112 return matrix;
113 }
114 //**********************************************************************************************
115
116 //**********************************************************************************************
124 inline const UpperMatrix<MT,SO,DF> generate( size_t n, size_t nonzeros ) const
125 {
128
129 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
130 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
131 }
132
133 UpperMatrix<MT,SO,DF> matrix( n );
134 randomize( matrix, nonzeros );
135
136 return matrix;
137 }
138 //**********************************************************************************************
139
140 //**********************************************************************************************
147 template< typename Arg > // Min/max argument type
148 inline const UpperMatrix<MT,SO,DF> generate( const Arg& min, const Arg& max ) const
149 {
151
152 UpperMatrix<MT,SO,DF> matrix;
153 randomize( matrix, min, max );
154 return matrix;
155 }
156 //**********************************************************************************************
157
158 //**********************************************************************************************
166 template< typename Arg > // Min/max argument type
167 inline const UpperMatrix<MT,SO,DF> generate( size_t n, const Arg& min, const Arg& max ) const
168 {
170
171 UpperMatrix<MT,SO,DF> matrix( n );
172 randomize( matrix, min, max );
173 return matrix;
174 }
175 //**********************************************************************************************
176
177 //**********************************************************************************************
187 template< typename Arg > // Min/max argument type
188 inline const UpperMatrix<MT,SO,DF>
189 generate( size_t n, size_t nonzeros, const Arg& min, const Arg& max ) const
190 {
193
194 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
195 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
196 }
197
198 UpperMatrix<MT,SO,DF> matrix( n );
199 randomize( matrix, nonzeros, min, max );
200
201 return matrix;
202 }
203 //**********************************************************************************************
204
205 //**********************************************************************************************
211 inline void randomize( UpperMatrix<MT,SO,DF>& matrix ) const
212 {
213 randomize( matrix, typename IsDenseMatrix<MT>::Type() );
214 }
215 //**********************************************************************************************
216
217 //**********************************************************************************************
225 inline void randomize( UpperMatrix<MT,false,DF>& matrix, size_t nonzeros ) const
226 {
228
229 using ET = ElementType_t<MT>;
230
231 const size_t n( matrix.rows() );
232
233 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
234 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
235 }
236
237 if( n == 0UL ) return;
238
239 matrix.reset();
240 matrix.reserve( nonzeros );
241
242 std::vector<size_t> dist( n );
243
244 for( size_t nz=0UL; nz<nonzeros; ) {
245 const size_t index = rand<size_t>( 0UL, n-1UL );
246 if( dist[index] == n - index ) continue;
247 ++dist[index];
248 ++nz;
249 }
250
251 for( size_t i=0UL; i<n; ++i ) {
252 const Indices<size_t> indices( i, n-1UL, dist[i] );
253 for( size_t j : indices ) {
254 matrix.append( i, j, rand<ET>() );
255 }
256 matrix.finalize( i );
257 }
258 }
259 //**********************************************************************************************
260
261 //**********************************************************************************************
269 inline void randomize( UpperMatrix<MT,true,DF>& matrix, size_t nonzeros ) const
270 {
272
273 using ET = ElementType_t<MT>;
274
275 const size_t n( matrix.rows() );
276
277 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
278 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
279 }
280
281 if( n == 0UL ) return;
282
283 matrix.reset();
284 matrix.reserve( nonzeros );
285
286 std::vector<size_t> dist( n );
287
288 for( size_t nz=0UL; nz<nonzeros; ) {
289 const size_t index = rand<size_t>( 0UL, n-1UL );
290 if( dist[index] == index+1UL ) continue;
291 ++dist[index];
292 ++nz;
293 }
294
295 for( size_t j=0UL; j<n; ++j ) {
296 const Indices<size_t> indices( 0UL, j, dist[j] );
297 for( size_t i : indices ) {
298 matrix.append( i, j, rand<ET>() );
299 }
300 matrix.finalize( j );
301 }
302 }
303 //**********************************************************************************************
304
305 //**********************************************************************************************
313 template< typename Arg > // Min/max argument type
314 inline void randomize( UpperMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max ) const
315 {
316 randomize( matrix, min, max, typename IsDenseMatrix<MT>::Type() );
317 }
318 //**********************************************************************************************
319
320 //**********************************************************************************************
330 template< typename Arg > // Min/max argument type
331 inline void randomize( UpperMatrix<MT,false,DF>& matrix,
332 size_t nonzeros, const Arg& min, const Arg& max ) const
333 {
335
336 using ET = ElementType_t<MT>;
337
338 const size_t n( matrix.rows() );
339
340 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
341 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
342 }
343
344 if( n == 0UL ) return;
345
346 matrix.reset();
347 matrix.reserve( nonzeros );
348
349 std::vector<size_t> dist( n );
350
351 for( size_t nz=0UL; nz<nonzeros; ) {
352 const size_t index = rand<size_t>( 0UL, n-1UL );
353 if( dist[index] == n - index ) continue;
354 ++dist[index];
355 ++nz;
356 }
357
358 for( size_t i=0UL; i<n; ++i ) {
359 const Indices<size_t> indices( i, n-1UL, dist[i] );
360 for( size_t j : indices ) {
361 matrix.append( i, j, rand<ET>( min, max ) );
362 }
363 matrix.finalize( i );
364 }
365 }
366 //**********************************************************************************************
367
368 //**********************************************************************************************
378 template< typename Arg > // Min/max argument type
379 inline void randomize( UpperMatrix<MT,true,DF>& matrix,
380 size_t nonzeros, const Arg& min, const Arg& max ) const
381 {
383
384 using ET = ElementType_t<MT>;
385
386 const size_t n( matrix.rows() );
387
388 if( nonzeros > UpperMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
389 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
390 }
391
392 if( n == 0UL ) return;
393
394 matrix.reset();
395 matrix.reserve( nonzeros );
396
397 std::vector<size_t> dist( n );
398
399 for( size_t nz=0UL; nz<nonzeros; ) {
400 const size_t index = rand<size_t>( 0UL, n-1UL );
401 if( dist[index] == index+1UL ) continue;
402 ++dist[index];
403 ++nz;
404 }
405
406 for( size_t j=0UL; j<n; ++j ) {
407 const Indices<size_t> indices( 0UL, j, dist[j] );
408 for( size_t i : indices ) {
409 matrix.append( i, j, rand<ET>( min, max ) );
410 }
411 matrix.finalize( j );
412 }
413 }
414 //**********************************************************************************************
415
416 private:
417 //**********************************************************************************************
423 inline void randomize( UpperMatrix<MT,SO,DF>& matrix, TrueType ) const
424 {
426
427 using ET = ElementType_t<MT>;
428
429 const size_t n( matrix.rows() );
430
431 for( size_t i=0UL; i<n; ++i ) {
432 for( size_t j=i; j<n; ++j ) {
433 matrix(i,j) = rand<ET>();
434 }
435 }
436 }
437 //**********************************************************************************************
438
439 //**********************************************************************************************
445 inline void randomize( UpperMatrix<MT,SO,DF>& matrix, FalseType ) const
446 {
448
449 const size_t n( matrix.rows() );
450
451 if( n == 0UL ) return;
452
453 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.3*n*n ) ) );
454
455 randomize( matrix, nonzeros );
456 }
457 //**********************************************************************************************
458
459 //**********************************************************************************************
467 template< typename Arg > // Min/max argument type
468 inline void randomize( UpperMatrix<MT,SO,DF>& matrix,
469 const Arg& min, const Arg& max, TrueType ) const
470 {
472
473 using ET = ElementType_t<MT>;
474
475 const size_t n( matrix.rows() );
476
477 for( size_t i=0UL; i<n; ++i ) {
478 for( size_t j=i; j<n; ++j ) {
479 matrix(i,j) = rand<ET>( min, max );
480 }
481 }
482 }
483 //**********************************************************************************************
484
485 //**********************************************************************************************
493 template< typename Arg > // Min/max argument type
494 inline void randomize( UpperMatrix<MT,SO,DF>& matrix,
495 const Arg& min, const Arg& max, FalseType ) const
496 {
498
499 const size_t n( matrix.rows() );
500
501 if( n == 0UL ) return;
502
503 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.3*n*n ) ) );
504
505 randomize( matrix, nonzeros, min, max );
506 }
507 //**********************************************************************************************
508};
510//*************************************************************************************************
511
512
513
514
515//=================================================================================================
516//
517// MAKE FUNCTIONS
518//
519//=================================================================================================
520
521//*************************************************************************************************
528template< typename MT // Type of the adapted matrix
529 , bool SO // Storage order of the adapted matrix
530 , bool DF > // Density flag
531void makeSymmetric( UpperMatrix<MT,SO,DF>& matrix )
532{
533 const size_t n( matrix.rows() );
534
535 reset( matrix );
536
537 for( size_t i=0UL; i<n; ++i ) {
538 matrix(i,i) = rand< ElementType_t<MT> >();
539 }
540
541 BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
542}
544//*************************************************************************************************
545
546
547//*************************************************************************************************
556template< typename MT // Type of the adapted matrix
557 , bool SO // Storage order of the adapted matrix
558 , bool DF // Density flag
559 , typename Arg > // Min/max argument type
560void makeSymmetric( UpperMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max )
561{
562 using Type = ElementType_t<MT>;
563
564 const size_t n( matrix.rows() );
565
566 reset( matrix );
567
568 for( size_t i=0UL; i<n; ++i ) {
569 matrix(i,i) = rand<Type>( min, max );
570 }
571
572 BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
573}
575//*************************************************************************************************
576
577
578//*************************************************************************************************
585template< typename MT // Type of the adapted matrix
586 , bool SO // Storage order of the adapted matrix
587 , bool DF > // Density flag
588void makeHermitian( UpperMatrix<MT,SO,DF>& matrix )
589{
590 using Type = UnderlyingBuiltin_t< ElementType_t<MT> >;
591
593
594 const size_t n( matrix.rows() );
595
596 reset( matrix );
597
598 for( size_t i=0UL; i<n; ++i ) {
599 matrix(i,i) = rand<Type>();
600 }
601
602 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
603}
605//*************************************************************************************************
606
607
608//*************************************************************************************************
617template< typename MT // Type of the adapted matrix
618 , bool SO // Storage order of the adapted matrix
619 , bool DF // Density flag
620 , typename Arg > // Min/max argument type
621void makeHermitian( UpperMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max )
622{
623 using Type = UnderlyingBuiltin_t< ElementType_t<MT> >;
624
626
627 const size_t n( matrix.rows() );
628
629 reset( matrix );
630
631 for( size_t i=0UL; i<n; ++i ) {
632 matrix(i,i) = rand<Type>( min, max );
633 }
634
635 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
636}
638//*************************************************************************************************
639
640
641//*************************************************************************************************
648template< typename MT // Type of the adapted matrix
649 , bool SO // Storage order of the adapted matrix
650 , bool DF > // Density flag
651void makePositiveDefinite( UpperMatrix<MT,SO,DF>& matrix )
652{
653 makeHermitian( matrix );
654}
656//*************************************************************************************************
657
658} // namespace blaze
659
660#endif
Header file for auxiliary alias declarations.
Header file for all basic DenseMatrix functionality.
Header file for the Indices class.
Header file for the IntegralConstant class template.
Header file for the IsDenseMatrix type trait.
Header file for the complete LowerMatrix implementation.
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 an upper 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) 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
bool isSymmetric(const DenseMatrix< MT, SO > &dm)
Checks if the given dense matrix is symmetric.
Definition: DenseMatrix.h:1456
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
constexpr void reset(Matrix< MT, SO > &matrix)
Resetting the given matrix.
Definition: Matrix.h:806
void randomize(T &&value)
Randomization of a given variable.
Definition: Random.h:626
#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.
Implementation of a random number generator.
Header file for basic type definitions.