Blaze 3.9
LowerMatrix.h
Go to the documentation of this file.
1//=================================================================================================
33//=================================================================================================
34
35#ifndef _BLAZE_MATH_LOWERMATRIX_H_
36#define _BLAZE_MATH_LOWERMATRIX_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< LowerMatrix<MT,SO,DF> >
83{
84 public:
85 //**********************************************************************************************
90 inline const LowerMatrix<MT,SO,DF> generate() const
91 {
93
94 LowerMatrix<MT,SO,DF> matrix;
95 randomize( matrix );
96 return matrix;
97 }
98 //**********************************************************************************************
99
100 //**********************************************************************************************
106 inline const LowerMatrix<MT,SO,DF> generate( size_t n ) const
107 {
109
110 LowerMatrix<MT,SO,DF> matrix( n );
111 randomize( matrix );
112 return matrix;
113 }
114 //**********************************************************************************************
115
116 //**********************************************************************************************
124 inline const LowerMatrix<MT,SO,DF> generate( size_t n, size_t nonzeros ) const
125 {
128
129 if( nonzeros > LowerMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
130 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
131 }
132
133 LowerMatrix<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 LowerMatrix<MT,SO,DF> generate( const Arg& min, const Arg& max ) const
149 {
151
152 LowerMatrix<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 LowerMatrix<MT,SO,DF>
168 generate( size_t n, const Arg& min, const Arg& max ) const
169 {
171
172 LowerMatrix<MT,SO,DF> 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 LowerMatrix<MT,SO,DF>
190 generate( size_t n, size_t nonzeros, const Arg& min, const Arg& max ) const
191 {
194
195 if( nonzeros > LowerMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
196 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
197 }
198
199 LowerMatrix<MT,SO,DF> matrix( n );
200 randomize( matrix, nonzeros, min, max );
201
202 return matrix;
203 }
204 //**********************************************************************************************
205
206 //**********************************************************************************************
212 inline void randomize( LowerMatrix<MT,SO,DF>& matrix ) const
213 {
214 randomize( matrix, typename IsDenseMatrix<MT>::Type() );
215 }
216 //**********************************************************************************************
217
218 //**********************************************************************************************
226 inline void randomize( LowerMatrix<MT,false,DF>& 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 > LowerMatrix<MT,SO,DF>::maxNonZeros( 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 );
242
243 std::vector<size_t> dist( n );
244
245 for( size_t nz=0UL; nz<nonzeros; ) {
246 const size_t index = rand<size_t>( 0UL, n-1UL );
247 if( dist[index] == index+1UL ) continue;
248 ++dist[index];
249 ++nz;
250 }
251
252 for( size_t i=0UL; i<n; ++i ) {
253 const Indices<size_t> indices( 0UL, i, dist[i] );
254 for( size_t j : indices ) {
255 matrix.append( i, j, rand<ET>() );
256 }
257 matrix.finalize( i );
258 }
259 }
260 //**********************************************************************************************
261
262 //**********************************************************************************************
270 inline void randomize( LowerMatrix<MT,true,DF>& matrix, size_t nonzeros ) const
271 {
273
274 using ET = ElementType_t<MT>;
275
276 const size_t n( matrix.rows() );
277
278 if( nonzeros > LowerMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
279 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
280 }
281
282 if( n == 0UL ) return;
283
284 matrix.reset();
285 matrix.reserve( nonzeros );
286
287 std::vector<size_t> dist( n );
288
289 for( size_t nz=0UL; nz<nonzeros; ) {
290 const size_t index = rand<size_t>( 0UL, n-1UL );
291 if( dist[index] == n - index ) continue;
292 ++dist[index];
293 ++nz;
294 }
295
296 for( size_t j=0UL; j<n; ++j ) {
297 const Indices<size_t> indices( j, n-1UL, dist[j] );
298 for( size_t i : indices ) {
299 matrix.append( i, j, rand<ET>() );
300 }
301 matrix.finalize( j );
302 }
303 }
304 //**********************************************************************************************
305
306 //**********************************************************************************************
314 template< typename Arg > // Min/max argument type
315 inline void randomize( LowerMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max ) const
316 {
317 randomize( matrix, min, max, typename IsDenseMatrix<MT>::Type() );
318 }
319 //**********************************************************************************************
320
321 //**********************************************************************************************
331 template< typename Arg > // Min/max argument type
332 inline void randomize( LowerMatrix<MT,false,DF>& matrix,
333 size_t nonzeros, const Arg& min, const Arg& max ) const
334 {
336
337 using ET = ElementType_t<MT>;
338
339 const size_t n( matrix.rows() );
340
341 if( nonzeros > LowerMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
342 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
343 }
344
345 if( n == 0UL ) return;
346
347 matrix.reset();
348 matrix.reserve( nonzeros );
349
350 std::vector<size_t> dist( n );
351
352 for( size_t nz=0UL; nz<nonzeros; ) {
353 const size_t index = rand<size_t>( 0UL, n-1UL );
354 if( dist[index] == index+1UL ) continue;
355 ++dist[index];
356 ++nz;
357 }
358
359 for( size_t i=0UL; i<n; ++i ) {
360 const Indices<size_t> indices( 0UL, i, dist[i] );
361 for( size_t j : indices ) {
362 matrix.append( i, j, rand<ET>( min, max ) );
363 }
364 matrix.finalize( i );
365 }
366 }
367 //**********************************************************************************************
368
369 //**********************************************************************************************
379 template< typename Arg > // Min/max argument type
380 inline void randomize( LowerMatrix<MT,true,DF>& matrix,
381 size_t nonzeros, const Arg& min, const Arg& max ) const
382 {
384
385 using ET = ElementType_t<MT>;
386
387 const size_t n( matrix.rows() );
388
389 if( nonzeros > LowerMatrix<MT,SO,DF>::maxNonZeros( n ) ) {
390 BLAZE_THROW_INVALID_ARGUMENT( "Invalid number of non-zero elements" );
391 }
392
393 if( n == 0UL ) return;
394
395 matrix.reset();
396 matrix.reserve( nonzeros );
397
398 std::vector<size_t> dist( n );
399
400 for( size_t nz=0UL; nz<nonzeros; ) {
401 const size_t index = rand<size_t>( 0UL, n-1UL );
402 if( dist[index] == n - index ) continue;
403 ++dist[index];
404 ++nz;
405 }
406
407 for( size_t j=0UL; j<n; ++j ) {
408 const Indices<size_t> indices( j, n-1UL, dist[j] );
409 for( size_t i : indices ) {
410 matrix.append( i, j, rand<ET>( min, max ) );
411 }
412 matrix.finalize( j );
413 }
414 }
415 //**********************************************************************************************
416
417 private:
418 //**********************************************************************************************
424 inline void randomize( LowerMatrix<MT,SO,DF>& matrix, TrueType ) const
425 {
427
428 using ET = ElementType_t<MT>;
429
430 const size_t n( matrix.rows() );
431
432 for( size_t i=0UL; i<n; ++i ) {
433 for( size_t j=0UL; j<=i; ++j ) {
434 matrix(i,j) = rand<ET>();
435 }
436 }
437 }
438 //**********************************************************************************************
439
440 //**********************************************************************************************
446 inline void randomize( LowerMatrix<MT,SO,DF>& matrix, FalseType ) const
447 {
449
450 const size_t n( matrix.rows() );
451
452 if( n == 0UL ) return;
453
454 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.3*n*n ) ) );
455
456 randomize( matrix, nonzeros );
457 }
458 //**********************************************************************************************
459
460 //**********************************************************************************************
468 template< typename Arg > // Min/max argument type
469 inline void randomize( LowerMatrix<MT,SO,DF>& matrix,
470 const Arg& min, const Arg& max, TrueType ) const
471 {
473
474 using ET = ElementType_t<MT>;
475
476 const size_t n( matrix.rows() );
477
478 for( size_t i=0UL; i<n; ++i ) {
479 for( size_t j=0UL; j<=i; ++j ) {
480 matrix(i,j) = rand<ET>( min, max );
481 }
482 }
483 }
484 //**********************************************************************************************
485
486 //**********************************************************************************************
494 template< typename Arg > // Min/max argument type
495 inline void randomize( LowerMatrix<MT,SO,DF>& matrix,
496 const Arg& min, const Arg& max, FalseType ) const
497 {
499
500 const size_t n( matrix.rows() );
501
502 if( n == 0UL ) return;
503
504 const size_t nonzeros( rand<size_t>( 1UL, std::ceil( 0.3*n*n ) ) );
505
506 randomize( matrix, nonzeros, min, max );
507 }
508 //**********************************************************************************************
509};
511//*************************************************************************************************
512
513
514
515
516//=================================================================================================
517//
518// MAKE FUNCTIONS
519//
520//=================================================================================================
521
522//*************************************************************************************************
529template< typename MT // Type of the adapted matrix
530 , bool SO // Storage order of the adapted matrix
531 , bool DF > // Density flag
532void makeSymmetric( LowerMatrix<MT,SO,DF>& matrix )
533{
534 const size_t n( matrix.rows() );
535
536 reset( matrix );
537
538 for( size_t i=0UL; i<n; ++i ) {
539 matrix(i,i) = rand< ElementType_t<MT> >();
540 }
541
542 BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
543}
545//*************************************************************************************************
546
547
548//*************************************************************************************************
557template< typename MT // Type of the adapted matrix
558 , bool SO // Storage order of the adapted matrix
559 , bool DF // Density flag
560 , typename Arg > // Min/max argument type
561void makeSymmetric( LowerMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max )
562{
563 using Type = ElementType_t<MT>;
564
565 const size_t n( matrix.rows() );
566
567 reset( matrix );
568
569 for( size_t i=0UL; i<n; ++i ) {
570 matrix(i,i) = rand<Type>( min, max );
571 }
572
573 BLAZE_INTERNAL_ASSERT( isSymmetric( matrix ), "Non-symmetric matrix detected" );
574}
576//*************************************************************************************************
577
578
579//*************************************************************************************************
586template< typename MT // Type of the adapted matrix
587 , bool SO // Storage order of the adapted matrix
588 , bool DF > // Density flag
589void makeHermitian( LowerMatrix<MT,SO,DF>& matrix )
590{
591 using Type = UnderlyingBuiltin_t< ElementType_t<MT> >;
592
594
595 const size_t n( matrix.rows() );
596
597 reset( matrix );
598
599 for( size_t i=0UL; i<n; ++i ) {
600 matrix(i,i) = rand<Type>();
601 }
602
603 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
604}
606//*************************************************************************************************
607
608
609//*************************************************************************************************
618template< typename MT // Type of the adapted matrix
619 , bool SO // Storage order of the adapted matrix
620 , bool DF // Density flag
621 , typename Arg > // Min/max argument type
622void makeHermitian( LowerMatrix<MT,SO,DF>& matrix, const Arg& min, const Arg& max )
623{
624 using Type = UnderlyingBuiltin_t< ElementType_t<MT> >;
625
627
628 const size_t n( matrix.rows() );
629
630 reset( matrix );
631
632 for( size_t i=0UL; i<n; ++i ) {
633 matrix(i,i) = rand<Type>( min, max );
634 }
635
636 BLAZE_INTERNAL_ASSERT( isHermitian( matrix ), "Non-Hermitian matrix detected" );
637}
639//*************************************************************************************************
640
641
642//*************************************************************************************************
649template< typename MT // Type of the adapted matrix
650 , bool SO // Storage order of the adapted matrix
651 , bool DF > // Density flag
652void makePositiveDefinite( LowerMatrix<MT,SO,DF>& matrix )
653{
654 makeHermitian( matrix );
655}
657//*************************************************************************************************
658
659} // namespace blaze
660
661#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.
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 complete UpperMatrix implementation.
Header file for the implementation of a diagonal matrix adaptor.
Header file for the implementation of a lower 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.