# This is a combination of 224 commits.

parent 6a18ebd545ce6eda1cb61fa5c3186865ffb722db
author Daniele G. Spampinato <[email protected]> 1713532895 +0200
committer Daniele G. Spampinato <[email protected]> 1713532921 +0200

# This is a combination of 178 commits.
# This is the 1st commit message:

Issue #303: framework for reference_dense added

This commit comes from a split of the original commit whose details are as follows:
SHA1 ID: a7c1db0eb45326184dadec39a6a7731eb16a39fd
Author: Albert-Jan Yzelman <[email protected]>  2022-01-14 15:52:26
Committer: Albert-Jan Yzelman <[email protected]>  2022-01-30 13:48:35

The reason for splitting is to separate the changes irrelevant to the branch 303 into a separate branch.
The split produced 3 separate commits.
This commit is number 1.
This commit in particular is relevant to 303 and will remain in it.

# This is the commit message #2:

Issue #303: unify script argument similar to --no-reference

# This is the commit message #3:

Issue #303: complete implementation of Vector< reference_dense >, modulo output iteration.

# This is the commit message #4:

Issue #303: modified vector implementation according to option 3 in handling dense container construction.

# This is the commit message #5:

Issue #303: vectors now have an iterator for the user to extract data from.

# This is the commit message #6:

Issue #303, fixed silly oops: the difference between two iterators should be an integer

# This is the commit message #7:

Issue #303: containers are now uninitialized on construction. Iteration over uninitialized containers returns nothing.

# This is the commit message #8:

Issue #303: configure and cmake now properly do make install with reference_dense

# This is the commit message #9:

Issue #303: cmake infra has updated-- make dense_reference additions compatible to it

# This is the commit message #10:

Drafting structured matrix interface using reference_dense backend

# This is the commit message #11:

Issue #303: Copy test mxm.cpp to dense_mxm.cpp as a starting point

# This is the commit message #12:

Issue #303: Add partial first version of Matrix interface for denseref

# This is the commit message #13:

Issue #303: WIP: dense mxm test

# This is the commit message #14:

Issue #303: Add dense mxm test to CMakeLists

# This is the commit message #15:

Issue #303: Separate reference and reference_dense includes

# This is the commit message #16:

Issue #303: Fix grb include name mismatch

# This is the commit message #17:

Issue #303: Add src files for reference dense backend

# This is the commit message #18:

Issue #303: Add project related option for dense backend

# This is the commit message #19:

Issue #303: Remove dense_reference-related code from reference cmake

# This is the commit message #20:

Keeping a clean separate denseref backend when building with cmake.

# This is the commit message #21:

Fix cmake option name

# This is the commit message #22:

Use a proper data container for a denseref matrix

# This is the commit message #23:

Change PRIVATE to INTERFACE for denseref backend

# This is the commit message #24:

fix cmake config for building backend

# This is the commit message #25:

Add friend functions for accessing matrix dimensions and raw data

# This is the commit message #26:

Add first version of blas3 mxm for reference dense backend

# This is the commit message #27:

Add temporary code for matrix data initialization

# This is the commit message #28:

Complete first working dense mxm unit test

# This is the commit message #29:

Added structured matrix test including draft static views and references

# This is the commit message #30:

Convert spaces to tabs

# This is the commit message #31:

Make use of operators and monoid provided to mxm

# This is the commit message #32:

Add friend getter/setter for field 'initialized' of denseref matrix

# This is the commit message #33:

Revert "Add temporary code for matrix data initialization"

This reverts commit 676cb4db7fe0ff54029c4b3b0ea31af60f8bc334.

# This is the commit message #34:

Add temporary way for allocating data of a denseref matrix

# This is the commit message #35:

Add functions to initialize matrix data from provided forward iterator

# This is the commit message #36:

Make getters const

# This is the commit message #37:

Exploit new features in mxm unit test

New features are building a matrix from an iterator and handling of
uninitialized matrices.

# This is the commit message #38:

added get/remove_ref type traits for StructuredMatrix types

# This is the commit message #39:

Temporary switch to pointers to target matrices to enable empty init.

# This is the commit message #40:

First draft of index mapping functions

# This is the commit message #41:

Cleaned up StructuredMatrix-related doxygen drafts

# This is the commit message #42:

Added to the spec draft

# This is the commit message #43:

Add alloc primitives for denseref

Currently they are a copy of the same primitives for reference backend.
Think whether it makes sense to group the alloc functions for these two
backends.

# This is the commit message #44:

Handle memory allocation failure. Implement destructor.

# This is the commit message #45:

WIP: Add data size calculator depending on Storage scheme and structure.

Also add a simple test to showcase the usage of this feature.

# This is the commit message #46:

Add buildMatrix for general full StructuredMatrix

# This is the commit message #47:

Add a check for size mismatch in buildMatrix

# This is the commit message #48:

Add mxm for dense Structured Matrices

# This is the commit message #49:

Introducing grb::get_view<> draft

# This is the commit message #50:

Add UpperTriangular StructuredMatrix (container and view)

# This is the commit message #51:

Implement polymorphic comparison of IMFs

# This is the commit message #52:

Add to structures a compatibility check from source structure and IMF functions

This is to be used when creating a view over a source matrix using
provided IMF functions. The compatibility of TargetStructure depends on
the SourceStructure and provided IMF functions. This is a runtime check.

# This is the commit message #53:

Add constructor taking two IMFs to reference SM with general structure

# This is the commit message #54:

Implement gather through view with provided TargetStructure and IMFs

# This is the commit message #55:

Add test for gathers on StructuredMatrix using views with IMFs

# This is the commit message #56:

Basic VectorView definition based on grb::vector container

# This is the commit message #57:

Draft of diagonal (vector) view of a matrix

# This is the commit message #58:

Initial get_view for VectorViews

# This is the commit message #59:

Add a stub for Storage Mapping Functions

# This is the commit message #60:

Define interface for dense StructuredMatrices and scalars

# This is the commit message #61:

Add test showcasing eWiseApply operations on Structured Matrices

# This is the commit message #62:

Minor changes for dense reference vector

# This is the commit message #63:

Add outer product of two vectors for reference_dense backend

# This is the commit message #64:

Expand the gather via IMF to general view of general Storage Matrix

# This is the commit message #65:

Add more Structures for StructuredMatrix

# This is the commit message #66:

Add dot product and norm2 in blas1 for dense reference backend

Also add a stub for a test showcasing their use

# This is the commit message #67:

WIP Constant matrix factory methods

# This is the commit message #68:

WIP Base classes for StructuredMatrix

# This is the commit message #69:

WIP Refactor StructuredMatrix using inheritance from Base classes

# This is the commit message #70:

Implement friend functions through base classes

This avoids the need to declare friend functions inside every
StructuredMatrix specialization.

# This is the commit message #71:

Add getter and setter for initialized field in StructuredMatrix

# This is the commit message #72:

Revert "WIP Constant matrix factory methods"

This reverts commit 742fe9f3ded2aaff5c389c6cdc6d6483081f756b.

# This is the commit message #73:

Add interface for creating constant matrices

This includes Identity matrix, Zero matrix and Givens rotations

# This is the commit message #74:

Add test showcasing the use of Identity and Zero matrices

# This is the commit message #75:

Revert "WIP: Add data size calculator depending on Storage scheme and structure."

This reverts commit 1f7b6c4f59c6473b7696a77ebb8b7d7b51358386.

# This is the commit message #76:

draft of selection and range views.

# This is the commit message #77:

Using default coords template param for VectorView.

# This is the commit message #78:

Added matrix-related doxygen draft. In need of compilation check (had to push to from machine without compiler support).

# This is the commit message #79:

Fixed interface compilation (need to double check buildMatrix friend functions though)

# This is the commit message #80:

Add doxygen for reference_dense vector

Also minor code cleanup.

# This is the commit message #81:

WIP blas1 interface for reference dense

# This is the commit message #82:

WIP: Add blas2 interface for reference_dense backend

# This is the commit message #83:

Fix compilation errors in denseref blas1

# This is the commit message #84:

Complete blas2 interface for denseref backend

# This is the commit message #85:

Fix compilation of vector

# This is the commit message #86:

Add possibility for VectorView to be a temporary

# This is the commit message #87:

Add VectorView tmp parameter to blas1 and blas2

# This is the commit message #88:

Add doxygen for VectorView

Provides detailed description of relevant concepts.

# This is the commit message #89:

Fix specification of tmp template parameters in blas1 and blas2

# This is the commit message #90:

Fix remaining tmp template errors in blas1

# This is the commit message #91:

Correct tmp template specializaiton in denseref VectorView

# This is the commit message #92:

Add missing tmp parameter to view source type

# This is the commit message #93:

Expand doxygen for denseref blas1

# This is the commit message #94:

Comment out performance semantics doxygen from blas1

It requires revision if still applies to the denseref backend

# This is the commit message #95:

Started working on doxygen description for level3 functions

# This is the commit message #96:

Added doxygen draft to blas3 operations

# This is the commit message #97:

compiler messages cleanup

# This is the commit message #98:

Add capacity to dense vector

# This is the commit message #99:

Expand doxygen for dense vector

# This is the commit message #100:

Inline functions to avoid "multiple definitions" linking error

# This is the commit message #101:

Remove tmp template parameter from VectorView interface.

This commit is the result of squashing the commits listed above and
cleanup of other locations where tmp parameter was introduced (doxygen,
blas3).

Revert "Add missing tmp parameter to view source type"
This reverts commit cd7099968ef5f0d72d9e2dcda00a93ee3ee1ad97.

Revert "Correct tmp template specializaiton in denseref VectorView"
This reverts commit bdc728e467c1fecb8ca9dd54f7b2aefc3af8d6b1.

Revert "Fix remaining tmp template errors in blas1"
This reverts commit 882183057fc0ec8e2578c553dbbc9f1b500bcaab.

Revert "Fix specification of tmp template parameters in blas1 and blas2"
This reverts commit d00b3edcba14275a86af6889f2d927d3daf73c88.

Revert "Add VectorView tmp parameter to blas1 and blas2"
This reverts commit d009aacdccb5d01f98c8bd32759c9c3321122a20.

Revert "Add possibility for VectorView to be a temporary"
This reverts commit fd5bf70cbceccf25595b5da1590a3f8b17f6c91a.

# This is the commit message #102:

Rename Identity view to Original

including related variables, typenames and doxygen fixes

# This is the commit message #103:

Rename StructuredMatrix base classes

This includes StructuredMatrixBase, StructuredMatrixContainer and
StructuredMatrixReference

# This is the commit message #104:

Remove tmp template parameter from StructuredMatrix interface

# This is the commit message #105:

Remove 'initialized' field from StructureMatrix classes

as it is already present in the underlying container grb::Matrix.
Expose that field through friend functions operaiting on
StructuredMatrix classes.

# This is the commit message #106:

Store reference instead of a pointer to a source matrix

# This is the commit message #107:

Adjust comments to reflect renaming Identity to Original view

# This is the commit message #108:

Remove 'initialized' field from VectorView

as it is already present in the underlying container grb::Vector.
Expose that field through friend functions operaiting on
VectorView.

# This is the commit message #109:

!19 Resolve internal issue #360: Make VectorView "compatible" with StructuredMatrix in terms of template parameters
This commit mainly adds Structure to the VectorView and reorders VectorView's template parameters to follow the same order as in StructuredMatrix. Doxygen is updated to reflect these changes.

In addition, there are the following changes:
* Remove VectorView's Tmp parameter from blas1 doxygen
* Use consistent ordering of template name parameters across blas 1/2/3. The order is Output, Mask, Inputs.
* Fix minor unrelated issues
* Generalize over storage scheme type for VectorView in level 1/2/3 routines

# This is the commit message #110:

!18 Resolve internal issue #357: Include capacity argument to VectorView and StructuredMatrix

# This is the commit message #111:

!20 Resolve internal issue #361: Introduce Scalar container

# This is the commit message #112:

!25 Resolve internal issue #366: Scalar version of operations at level 1

This commit introduces the following changes:
* Update ALP tests to use ALP Scalar
* Update level 1 functions to operate on ALP scalars

# This is the commit message #113:

!30 Resolve internal issue #377: Do not expose coordinates template argument in alp vectors

# This is the commit message #114:

!24 Resolves internal issue #371: Introducing Union of Bands Structure
* Using past-the-end notation for static grb::Interval
* Indent comments from MR !24
* simple example of usage in test
* Initial draft of Band structure

# This is the commit message #115:

!26 Resolve internal issue #364: Resize operations for reference dense containers

* Add resize functions for dense ALP containers
* Resize for dense containers succeeds if new_cap <= current_cap

# This is the commit message #116:

!29 Resolve internal issue #362: Implement bracket operator for accessing single element

* Add bracket operator for ALP scalars.
* Add bracket operator to VectorView

# This is the commit message #117:

!31 Resolves internal issue #376: All views applied via get_view
* Dropped leftover coord arguments in VectorView
* Integrating comments from MR !31
* Doxygen cleanup
* Started using get_view for all views in matrix.hpp

# This is the commit message #118:

!23 Resolve internal issue #370: Introduce sort operation on vectors
* Move sort to level-1 functions
* clean buffer and move to C++ comparators
* define sort interface

# This is the commit message #119:

!32 Resolve internal issue #378: Convert grb storage from namespace to enum

This enum now has two values - Dense and Sparse.
The actual implementation of the storage scheme
will be handled by the Storage Mapping Functions.

# This is the commit message #120:

!35 Resolve internal issue #367: Separate alp and grb

This commit is part of two-commit sequence that resolves internal issue #367.
The commits are separated to allow for easier tracking of the moved files' history.

This includes:
 - Add alp cmake infrastructure
 - Add include/alp.hpp
 - Add necessary cpp files
 - Introduce include/alp and src/alp subfolders
 - Move include/graphblas/denseref to include/alp/reference
 - Move ALP-related root hpp files from include/graphblas to include/alp
 - Copy needed graphblas files to alp subfolders
   - some include root files
   - some include/base files
   - some include/utils files
 - Remove reference_dense backend implementation from graphblas as it is now part of ALP

# This is the commit message #121:

!35 Resolve internal issue #367: Refactor grb to alp

This commit is part of two-commit sequence that resolves internal issue #367.
The commits are separated to allow for easier tracking of the moved files' history.

This includes the following changes in relevant files
(include, src and tests):
 - Rename namespace from grb to alp
 - Rename backend from reference_dense to reference
 - Rename include guards
 - Remove graphblas specific backends from alp
 - Other minor related changes
 - Fix previously undiscovered compilation issues
 - Unrelated: Remove nrows/ncols from blas2. Move its doxygen to matrix.hpp

In addition, refactor naming of ALP containers
 - Rename Matrix -> internal::Matrix
 - Rename StructuredMatrix -> Matrix
 - Rename Vector -> internal::Vector
 - Rename VectorView -> Vector
 - Modifications incorrporating the above-mentioned renames
 - Minor coding style fixes where appropriate (mainly indentations)

# This is the commit message #122:

Add doxygen configuration file for public ALP API

# This is the commit message #123:

Add ALP doxygen main page. Also some minor code style fixes.

# This is the commit message #124:

!33 Resolve internal issue #379: Implement use cases with current ALP interface
* Add use cases for ALP interface
    - gemm
    - cholesky
    - householder tridiagonalization
    - conjugate gradient with sparse containers
* Cholesky not using diagonal view for alpha
* Other changes required to implement the algorithms:
    - Added trapezoidal structure
    - Added get_view with perm vector.
    - Moved constant structured matrices functions into the appropriate namespace
    - Add foldr( Scalar, Matrix, monoid )
    - Add row / column views on matrices

# This is the commit message #125:

!44 Resolve internal issue #395: Lazy execution of functions resulting in low-rank containers

This commit adds support for low-rank containers resulting from certain operations like outer product. These containers are internally implemented as views over lambda functions. Main goal of this design is to allow lazy calculation of resulting matrix's elements and avoid allocating the container. The views over lambda functions are internal functionality (not exposed to the user).

There are two new outer functions added:
 - one that computes x * y^T and returns a matrix with a general structure,
 - one that computes x * x^T and returns a matrix with a symmetric structure.

In addition, this commit
- Updates cholesky to use new interface
- Improves the existing unit test for outer product on dense containers

# This is the commit message #126:

!9 Resolve internal issue #410: Rework dense matrix container

This commit introduces the following changes:
- Implement Index Mapping Functions (IMFs) using static polymorphism.
- Introduce Storage Mapping Functions (SMF) that map coordinates
from logical space into the physical space (actual container in
memory). SMFs are modeled as quadratic bivariate functions, where
variables are matrix coordinates in logical space.
- Introduce Acess Mapping Function (AMF) which encapsulates two IMFs
(one per dimension) and SMF into a single function. Enable static
fusion of specific families of IMFs with the SMF.
- Refactor Matrix container to adopt above-mentioned changes. These
changes enable simplification of the matrix class hierarchy and reduce
the code size.

In addition, this commits makes the following minor unrelated changes:
- Improve interface and implementation of nrows, ncols and dims
- Expose the type of the base matrix
- Temporarily disable some unit tests until rework is completed
- Remove trailing whitespaces
- Adapt level2/3 and IO functions signatures to the new matrix signature

# This is the commit message #127:

Add a getter for a reference to an element in matrix (#16)


This introduces an internal interface for accessing elements of a matrix.
It is designed to be used by operations and not exposed to the user.

There are two functions, getStorageIndex and getCoords, that enable mapping from logical to physical iteration spaces and vice versa, respectively.
Finally, the third function, access, takes storage index (physical iteration space) as input and returns a reference or a value of the corresponding element, depending if a Matrix is a container matrix or a functor-like object.

# This is the commit message #128:

Resolve internal issue 434: AMF to provide physical container size (#17)


This commit introduces a new feature to AMF - providing storage size of the corresponding container,
which is then used to allocate the physical container. Currently, the storage size is provided to the AMF
by the matrix class as this is the only place where all relevant factors (e.g., structure, logical dimensions)
are known.

The changes relevant to the above are:
* Expose the sizes of container's logical and physical layouts
* Improve the design of Storage mapping functions
* Encapsulate AMF creation in a factory class

This commit also makes the following changes:
* Fix visibility of AMF's class members
* Improve naming of certain functions/types
* Improve doxygen and code style of storage related code.
* Rename smf.hpp into storage.hpp
* Rename namespace smf to storage

Unrelated fixes:
* Fix the type of IMF composition when creating views over matrices

# This is the commit message #129:

Resolve internal issue 445: update dense reference vector implementation (#20)


Implement dense ALP Vector as N x 1 dense ALP Matrix

The interface offers Vector as a separate class,
but behind the scenes, vector inherits from Matrix
and makes use of all its existing mechanisms.

Some relevant changes:
* Allow matrix classes to be inherited from
* Add IMF to vector's template parameter list
* Implement vector as a N x 1 specialization of dense ALP matrix
* Also, adapt reference::matrix to vector signature change.
* Adapt vector interface in operations according to new changes

In addition, this commit also introduces the following unrelated changes:
* Fix indentation in reference/vector.hpp
* Improve coding style (formatting) of template arguments (aligning, ordering, etc)
# This is the commit message #130:

412 dense alp mxm (#21)

* Reorganized scalar so to have backend-specific implementations. Added reference implementation.
* added basic, union-band mxm for union-of-band input and possibly symmetric output (symmetric input not yet supported)
# This is the commit message #131:

Resolve internal issue 448: internal set getinitialized to be invoked via matrixbase (#22)

* Implement get/setInitialized call through MatrixBase
* Add a const qualifier for boolean argument of existing setInitialized signatures.
* Add 'initialized' field to MatrixFunctor
* Add get/setInitialized to MatrixFunctor
* Annotate that initialized field in MatrixFunctor is a temporary solution
# This is the commit message #132:

455 revise matrix access call stack for dense alp (#23)

* Static cast to a reference instead of a type
* Fix code style in IMF classes
* Remove unnecessary const qualifier
# This is the commit message #133:

Rename internal intermediate classes to avoid confusion with existing concepts (#25)

MatrixContainer -> StorageBasedMatrix
MatrixFunctor -> FunctorBasedMatrix
# This is the commit message #134:

453 dense alp mxm symmetry (#26)

* Taking symmetry into account in mxm

* Adding numerical tests for structured mxm

* First numerical tests on ge, tr, sy output

* Changes according to review comments (0/1 integer variables to bool)
# This is the commit message #135:

458 container level traits (#24)

Add type traits to:
* check if a type is ALP container
* inspect structure and view types of ALP containers
* check if provided view type is associated with an ALP container type that requires memory allocation
* check if provided view type is a view over a storage-based ALP container
* check if provided view type is a view over a functor-based ALP container

Expose internally the type trait for checking if a type is an internal container (matrix/vector).

Add a unit test to showcase the use of the added type traits.

Unrelated changes:
* Fix the type of transposed view for a square matrix

# This is the commit message #136:

454 generalize matrix over container functor (#27)

* Change signature of MatrixFunctor class

Avoid depending on std::function and
pass lambda function type as the class template parameter.
Also, pass store function as an object instead of as a reference.
Pass lambda function by value in the constructor.

* Declare imfs and lambda as constant class members in MatrixFunctor

* Order MatrixFunctor constructor parameters consistently with MatrixContainer

* Expose lambda type in the functor matrix

It is more readable than View::applied_to

* Create Identity matrix as a MatrixFunctor

* Minor fixes

* Enable compilation of constant matrix unit test

* Adapt the constant matrix unit test to the current repo state

* Make MatrixBase a friend of MatrixFunctor so it can call its methods

* Restore identity matrix implementatin to use std::function instead of pure lambda

* Remove obsolete code and comments

* Generalize MatrixFunctor for constant matrices

* Adapt Zero matrix factory method to the current interface

* Enable set/getInitialized in reference/matrix.hpp only for ALP matrices

* Make Matrix conditionally inherit from Container of Functor

depending on whether it is a view over a container or a lambda function.

* Adapt Functor Matrix constructor calls to the new interface

* Fix the implementation of getInitialized for ALP vector

* Expose certain type definitions within ALP vector

* Fix rebase mistake

* Use view type traits to differentiate between different internal ALP matrix specializations

* Properly name the static property indicating whether a container allocates memory

Also use new type traits to calculate its value.

* Add type traits to check the underlying nature of an ALP container

Check if a type is storage-based ALP container, and
Check if a type is functor-based ALP container.

* Add const reference getter to a functor object of a functor-based matrix

* Improve the SFINAE conditions for storage-based matrix containers

* Add matrix constructors for functor and view based matrices.

Also update the SFINAE of functor non-view based constructor to use new type traits

* Update doxygen to reflect the SINAFE mechanism changes

* Add missing parameters to getStorageIndex in reference vector

This should be revised for parallel backends.

* Re-enable dense outer unit test

* Update other matrix specializations with the new constructor selection mechanism.

Also add all missing functor-based constructors to these specializations.

* Remove obsolete speciaization for constant matrices

The proper implementation for constant matrices is pending.
Disable constant matrices check until the constant matrix specializations are implemented.

* Code review fix: add Vector functor constructors

Also use new type traits for compile-time enabling of constructors.
Reorder existing constructors and update doxygen according to the order in Matrix class.

* Code review fix: Correct the return type of getFunctor forward declaration

* Remove unnecessary reference to size_t lambda function parameters

* Fix dims relatd issues in matrix base classes

* Correct the base type property of matrix

* Correct type in static_cast in get/setInitialized for vectors

* Add a temporary version of buildVector

* Update dense_outer test

* Change the way lambda function is used in functor-based matrices

* Use consistent SFINAE enable_ifs

* Fix the base type of various matrix specializations

* Add correctness check in the dense outer product unit test

Utility functions borrowed and adapted from dense mxm unit test.
This test only tests functor-based results. The storage-based result
cannot be tested because the corresponding blas function is not yet
implemented.
# This is the commit message #137:

450 initialized status of matrixfunctor (#29)

* The init status of functor-based containers runtime-depends on input containers
* Rename data-related lambda function to avoid confusion with init lambda
* Add functor init status check to outer unit test
# This is the commit message #138:

470 implement the change of amf on transposition (#30)

* Simplify template specification in some get_view functions
* Save templated values for static coefficients as class members
* Add static mapping of view type to the enumerator value
* Add a special view enum value for view types not exposed to users
* Implement type traits to transform polynomial type according to the provided view
* Add methods to create AMF by applying views onto provided AMF
* Add matrix and vector constructors for creating views with provided AMF
* Mapping polynomial type depends on the view type of the matrix
* Add unit test for views
* Temporarily disable test for diagonal views
# This is the commit message #139:

478 amffactory and matrix views (#31)

* Clean up AMF specialization for static IMFs

Store only Id IMFs to keep the size.
Original IMFs are already fused into the mapping polynomial.
Remove obsolete fields and make fused function static.

* Correct the visibility of AMF in storage based matrix

* Group the mechanisms for creating AMF objects. Encapsulate internals inside a single AMF type in StorageBasedMatrix.

Offer the following three ways to instantiate AMFs:
* AMF belonging to an original container (that requires allocation)
* AMF as a result of a gather view using row and column IMFs
* AMF as a result of a static view (e.g. transpose, diagonal)

Other related changes:
* Implement fusion of static IMFs into the storage mapping polynomial.
* Add a type trait to extract amf_type from any type (to support invocation for original containers)

Change the definition of StorageBasedMatrix:
* Encapsulate IMFR, IMFC and polynomial in a single AMF type
* Add type traits to determine the AMF type for a given matrix static properties.

* Explicitely set matrix init status to true in the test so that accesses can succeed

* Temporarily disable outer test until views on functor-based containers are functional

* Remove the unneeded AMF specialization for Strided IMFs

Since Strided IMFs can always be fused into the polynomial, this specialization will never be instantiated.

* Specialize getStorageIndex for cases when at least one of the IMFs is an Id

* Use structs for AMF factories

* Simplify the logic to determine Matrix AMF type

* Remove obsolete type traits

* Improve the handling of AMF objects

Instantiate original AMF only through factory
and disable making AMF objects outside of AMF factory
Disallow AMF copies and enforce move semantics

* Fix wrong order of parameters in composed IMF constructor

* Use more precise naming for AMFFactory classes

Also disable instantiation of these classes

* Avoid unnecessary computation of b and s parameters when composing with an ID IMF

* Statically check that transpose and diagonal view come with ID IMFs

* Improve naming of IMFs in composition

Also solve inconsistent ordering of imf-type function parameters

* Delete default constructor instead of making it private
# This is the commit message #140:

478 vector rework and vector views (#43)

This commit introduces the following major changes:
* The way ALP vector uses an ALP matrix is modified to allow for better support of views that result in vectors.
* Implement views on matrices and vectors that result in vectors:

The detailed list of changes is as follows:
* Use two IMFs in ALP Vector to allow for a better 1D-2D mapping
* Implement row/column/diagonal views over a matrix
* Update the naming of AMF factory classes
* Add a method to create AMFs for original matrix when used as a vector
* Move definition of FromPolynomial after Compose as the first uses the second
* Fix the order of IMFs when creating a IMF composition
* Simplify IMF composition type trait specializations
* Add Zero IMF and use it as a default column IMF for vector
* Add type trait to determine storage mapping polynomial based on matrix static properties
* Make IMF constructors explicit
* Use poly calcuation in the determine_amf_type
* Provide IMF do AMF factory instead of size to allow for constructing based on Id and Zero IMFs
* Add mapping polynomial specific to Vectors
* Implement gather views over vectors
* Expand testing of vector views
* Add backend as an input parameter to polynomial selection logic
* Access AMF-internal types directly and do not expose them through matrix
* Re-enable previously disabled lines in unit test for type traits
# This is the commit message #141:

Update Functor-based matrices according to the recent AMF changes (#51)

Also, re-enable and improve outer unit test
# This is the commit message #142:

493 simplify the design of alp containers (#53)

The detailed list of changes is as follows:
* Centralize Matrix base class selection logic
* Expose amf type just in storage based matrix
* Call the base constructor via the defined type
* Add missing constructors for Matrix view with provided AMF
* Improve SFINAE logic for Matrix constructors
* Improve SFINAE logic for vector constructors
* Remove amf_type declaration from vector as it is already defined in its base class
# This is the commit message #143:

497 fix propagation of imf types of views over a matrix vector (#54)

* Added Band Matrix type.

* Fixed propagation of IMF types via view_type

* Fix indentation
# This is the commit message #144:

499 Clean up get_view interface (#56)

This commit makes the following contributions:
* Move all AMFFactory invokations to Matrix constructors
* Add type trait that defines new container types with a possibility to modify a single template parameter.
* Add vector view type used for gather view on a matrix when one dimension is of size 1 (returning a vector)
* Remove redundant get_view implementations
* Add a gather view structure and enum value

Minor fixes:
* Rename target to source when refering to the source container of a view
* Clean up doxygen for get_view functions and relevant container constructors
* Also remove the now unnecessary vector view
* Remove deprecated code
* Mark matrix/vector constructor taking AMF as deprecated
# This is the commit message #145:

504 vector views on a matrix should result in a canonical vector representation (#57)

All vectors are now represented as column matrices.
The list of relevant changes:
* Fix remaining inconsistencies between original and gather views
* Reimplement diagonal view to reshape the polynomial by combining j coefficients into i coefficients
* Add constant IMF
* Implement column view over a matrix as a special case of gather view with a constant IMF
* Adapt vector-matrix interaction on access to follow the other relevant changes within this branch
* Reimplement row view over a matrix as a column view over a transposed source
* Fix return type of the simple gather view
* Improve constant IMF by allowing different domain and range
* Gather view over a vector always has ID as the column IMF
* Fix calculation of poly coefficients in AMF Reshape factories
* Reimplement diagonal view to be a two step process in case of non-square matrices
* Remove the type of diagonal view from general matrix definition
# This is the commit message #146:

512 implement matrix view over a vector (#60)

Taking a matrix view over a vector produces a column matrix of size Mx1, where M is the vector length. In other words, the vector is observed as a column matrix. As a result, vectors may be used as matrices in operations that operate on matrices.
The list of changes this commits introduces is the following:
* Add types to support Matrix view over a vector
* Add AMF reshape specialization for creating a matrix view
* Add get_view that implements a matrix view over a vector
* Extend unit test with the new view type
# This is the commit message #147:

type_trait draft for checking that a band_interval type contains or not another one.

Also fixing the following assumption that intervals should be defined compact: meaning two interval have to be at least at distance one from each others.

# This is the commit message #148:

Fixed typo

# This is the commit message #149:

print bands

# This is the commit message #150:

foldl( vector, scalar, op )

# This is the commit message #151:

temp fixes

# This is the commit message #152:

foldl(vec,scalar) and foldl(scalar,vec) with unit tests

-revert dense_dot_norm2
-foldr calls foldl internally

# This is the commit message #153:

revert non-fold related changes to 303

# This is the commit message #154:

missing unit test file

# This is the commit message #155:

review 01

# This is the commit message #156:

review 2

# This is the commit message #157:

494 implement set for alp containers scalar to matrix (#58)

Introduces the following changes:
* Implement set operation for initializing all elements of a matrix with a given scalar
* Add a test for set with ALP containers
* Add set( matrix, matrix ) implementation
* Fix the support for symmetry
* Limit iteration space over row coordinate
* Add band support to set(matrix, scalar)
* Fix wrong upper j limit and add declare limit variables as const
* Remove masked set variants
* Implement set on vectors
* Add missing foldl( vector, vector )

In addition to these, adds the following minor fixes in containers
  * Fix cv-qualifier for getInitialized( Scalar )
  * Add a method to get storage dimensions of a storage-based matrix
  * Add const reference value getter to scalar
  * Do not check for init-status when accessing internal::vector with non-const operator[]
  * Other minor fixes
# This is the commit message #158:

516 implement eWiseLambda (#63)

* Fix includes in reference/blas1.hpp
* Implement eWiseLambda for vectors and matrices
Vector variant creates a matrix view over a vector, wrapper for a lambda function that takes two coordinates instead of one, and finally calls matrix eWiseLambda variant with these parameters.
* Add unit test for eWiseLambda with simple numerical checks
# This is the commit message #159:

505 dense dot norm2 (#66)

* dot implemented via apply and foldr

rebased on foldr
only basic unit test is added

* uze setInitialized for z

* unit tests done

TO DO after rebase:
				enable set()
				enable foldr()
    enable dot()

* dot() check if scalar is initialized

* rebase and cleanup

* one more temp disable

* Fix support for functors in ALP vector

* Enable commented-out code in dot product and temporarily disable test 2 in the unit test

* remove redundant captures in lambda
enable all tests

* cleanup unit test

* review 01

Co-authored-by: Denis Jelovina <[email protected]>
Co-authored-by: Vladimir Dimic <[email protected]>
# This is the commit message #160:

480 storage scheme selection and instantiation (#67)

* Add formulation for packed format and clean up the doxygen for storage formats

* Implement factory classes for storage polynomials

Update ALP matrix with new Poly factory interface

* Add mechanism to select polynomial type based on static matrix properties

* Change AMF Factory FromPolynomial to take PolyFactory as template parameter

This encapsulates polynomial creation logic within the AMF factory and simplifies the invocation of AMF Factory Create function

* Change the type of dynamic poly coefficients to a signed integer type to allow negative coefficients

* Move storage dimension calculation to poly factories

* Add a definition of abstract poly factory.

This serves only as documentation.

* Reimplement Packed factories

* Add stub for banded factory

* Use enums as template parameters instead of bools for better readability

* Add unit test for storage polynomials
# This is the commit message #161:

Hotfix: bound calculations and output init status in set, return value in buildMatrix

# This is the commit message #162:

Hotfix: Add ALP matrix specialization for LowerTriangular structure

# This is the commit message #163:

Hotfix: Fix the structure of column/row get_view return type

# This is the commit message #164:

521 clean up utility functions for printing alp containers (#68)

* Add functionalities to check whether a given coordinate falls inside a non-zero area defined by a union of bands
* Add utility functions for printing ALP containers
* Use new print functions in unit tests
# This is the commit message #165:

523 dense file iterator (#70)

Implementation of file reader and iterator for use in buildMatrix() and similar.

This implementation reads only  dense matrices saved in Matrix Market Array Format.

Currently only forward access iterator is supported. 
Matrix Market Array Format stores elements in column-wise order and this iterator is returning elements in the same order.

For matrices with symmetry, Symmetric and  Hermitian matrices are supported, 
but iterator reads and returns only lower triangular part of matrix.
It should be taken care, within ALP code, to build the proper matrix, from this lower triangular part, or implement the corresponding logic.

Buffer for matrix elements has been removed, and should be implemented in the future.

It is suggested that this utility should be merged with analogous reader for sparse matrices, implemented in GRB namespace. 


Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #166:

Resolve internal issue #513: iterating within bands in operations (#69)

* Add functions to calculate row- and column-coordinate limits. This is used by operations to calculate iteration space for each band.
* Use new utility functions in operations. Also minor code style and documentation fixes in the modified functions.
* Improve naming of the functions that return band limits
* Improve the logic for iteration space calculation
# This is the commit message #167:

Resolve internal issue #519: Implement eWiseApply (#71)

* Implement eWiseApply for dense containers

Adds a version that operates on a single band.
Use ALP scalars instead of plain types as inputs to eWiseApply.
Remove allow_void in eWiseApply.
Improve readability and other minor code style fixes.

* Add unit test for eWiseApply operations and remove the old one
# This is the commit message #168:

520 dense cholesky integration (#72)

run like this:
make test_alp_cholesky_ndebug_alp_reference &&  tests/unit/alp_cholesky_ndebug_alp_reference  -fname datasets/spd5.mtx
or 
make test_alp_cholesky_ndebug_alp_reference &&  tests/unit/alp_cholesky_ndebug_alp_reference  -n 100

The first version reads  dense symmetric real matrix in Matrix Market  Array format,  the second version  generated a random SPD matrix internally. Once input matrix is present the code calculates Cholesky decomposition and checks the result. 


* views fix

* test cholesky

* First numerically correct version

* Add missing structure type traits

* Use utils print functions

* Add preliminary checks for cholesky result

* Fix a bug in print_matrix

* Fix band upper limits as the end is not included

TODO: Check other structures

* file reader enabled

currently works on general input matrices:
make test_alp_cholesky_ndebug_alp_reference &&  tests/unit/alp_cholesky_ndebug_alp_reference  -fname ../datasets/spd5_convertotodense.mtx

But it  should also work for symmetric input matrices
make test_alp_cholesky_ndebug_alp_reference &&  tests/unit/alp_cholesky_ndebug_alp_reference  -fname ../datasets/spd5_convertotodense_symmetric.mtx

* clean prints and DEBUG

* Unit test with check or the solution by
calculating the Frobenius norm of H - L L^t
where H is the input matrix and L is Cholesky solution

H is expcted to be (currently real) symmetric matrix in
Matrix Market Array format

The code is currently generating Upper instead of lower triangular matrix,
which is solved by using the transpose view.

* remove temporary parser in implement it in a separate branch

* review and internal matrix generation

* remove temp matrix files

* spaces

Co-authored-by: Denis Jelovina <[email protected]>
Co-authored-by: Vladimir Dimic <[email protected]>
# This is the commit message #169:

Resolve internal issue #408: Remove deprecated functions (#73)

* Remove deprecated eWiseMulAdd and eWiseAdd
* Replace delegating calls to deprecated functions with temporary return PANIC
# This is the commit message #170:

528 dense foldr explicit implementation (#77)

* foldr and some unit tests

* vector->scalar unit test

* foldl/r unit tests with views

* foldl/r vector->vector versions and unit tests

* review

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #171:

Resolve internal issue #526: Implement eWiseMul for matrices with matching structures (#76)

* Remove masked variants of eWiseMul as they are not supported in the dense backend
* Implement eWiseMul for matrices with matching structures
* Add test for eWiseMul on matrices
# This is the commit message #172:

Resolve internal issue #494: Implement set for scalars (#80)

* Add set for scalar

* Add the new set to the unit test
# This is the commit message #173:

Resolve internal issue #531: Implement foldl and foldr on matrices (#78)

* Add internal fold variants ( matrix <- matrix|scalar )
* Add foldl and foldr for ewise folding a scalar or a matrix into a matrix
* Expand fold unit test to include new matrix-based folds
* Add operator variants for fold on the user-level
* Add a TODO comment to check non-zero structure compatibility
# This is the commit message #174:

Fixes wrong symmetry constraints in mxm_band_generic (#81)


# This is the commit message #175:

Resolve internal issue #530: Integrate gemm-like example (#82)

This introduces a gemm-like algorithm implemented using ALP primitives and a test showcasing its use.

Some changes compared to the previous WIP version are:
* Test all combinations of transposed/non-transposed A and B
* Fix const qualifiers of input parameters in gemm_like_example.
* The test uses randomly generated matrices as inputs.
* The numerical correctness check verifies that the modified values are correct as well as that values that should not be modified are indeed not modified.
* Add dimension missmatch check in the gemm-like algorithm
* Fix view size when using stride > 1
* Use more general view parameters (stride and starting elements)
* Add comments and debug prints
# This is the commit message #176:

Resolve internal issue #541: Add ALP matrix specialization with an orthogonal structure (#84)

* Add orthogonal matrix. This includes specialization of alp::matrix and polynomial selection type traits
* Add orthogonal matrix to unit test
# This is the commit message #177:

Resolve internal issue #494: Implement matrix set via fold (#83)

* Implement set on matrices using foldl on matrices
* Remove obsolete internal set functions
* Fix code style issues of modified set operations
* Add TODO comment about structure match checking
# This is the commit message #178:

Resolve internal issue #542: Add ALP matrix specialization with symmetric tridiagonal structure (#86)

* Fix upper band limit for tridiagonal structures
* These structures infer properties from a specific band structure rather than open interval
* Break long lines
* Add matrix specialization for symmetric tridiagonal structure
* Extend the unit test with symmetric tridiagonal matrix and querying symmetry

# This is the commit message #179:

Fix after rebase onto develop: Scalar is not part of graphblas

# This is the commit message #180:

Fix after rebase onto develop: Remove duplicate functions

It seems that these were moved to graphblas/base

# This is the commit message #181:

Fix after rebase onto develop: Remove unneded function

# This is the commit message #182:

Resolve internal issue #547: Add hermitian structure and the corresponding ALP matrix specialization (#87)

* Add Hermitian structure
* Add a determine_poly_factory for hermitian structure
* Add matrix specialization for hermitian structure
* Expand unit test with a hermitian matrix
* Declare hermitian matrix with a complex data type
# This is the commit message #183:

Resolve internal issue #548: Add hermitian tridiagonal structure and the corresponding ALP matrix specialization (#88)

* Add hermitian tridiagonal structure
* Add specializations for storage polynomial selector and ALP matrix for Hermitian Tridiagonal structure
* Extend unit test with a hermitian tridiagonal matrix. Also extend ask_questions with a check if a matrix is tridiagonal
* Remove redundant structure from inferred structures
# This is the commit message #184:

Fix the use of set on scalars (#90)

Remove the temporary workaround.
# This is the commit message #185:

Resolve internal issue #549: Add conjugate operation for ALP matrix (#89)

* Add conjugate operation with specializations for square and non-square matrices
* Add unit test for conjugate. It checks the correctness of conjugate by comparing the result with the tranpose of the original matrix. Tests complex and real matrices.
* Add an explanation for a temporary use of square instead of symmetric structure for real matrices

# This is the commit message #186:

Resolve internal issue #553: Allow for a tolerance in error checking in ALP gemm unit test (#91)


# This is the commit message #187:

551 fix some script (#92)

* split old and new smoke tests
* move alp_cholesky to smoke and enable smoke test
* move alp_gemm to smoke and enable smoke test

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #188:

532 integrate householder with dense alp containers (#93)

* compile and test with
     -make test_alp_zhetrd_alp_reference && tests/smoke/alp_zhetrd_alp_reference
     -make test_alp_zhetrd_complex_alp_reference && tests/smoke/alp_zhetrd_complex_alp_reference
* numerical correctness checked by Frabenus norm H - QTQt
* additional check 
   -matrix Q check, if it is orthognal
   -matrix QTQt check, if it is a correct solution

Non-householder related fixes
* print matrices in python like style
* print matrices with max precision 1.e-10
* initialize identity using view::diagonal
* fix missing dot() and norm2() scalar conversion
* outer enable complex conjugate when complex numbers used
* conjugate function on matrices and vectors added
* fix the view-related type traits for Hermitian structure
* square view fix
* move alp_zhetrd to smoke tests
* add checks of alp_zhetrd in smoketests.sh (complex and real version)
* view() is General by default

issues
1) Hermitian currently using Square storage
2) TridiagonalHermitian  using Square storage
3) symmetric view of general, done using eWiseLambda
4) fold on non-matching-structures  Matrices does not work
5) transpose on Symmetric does not work when reused on Square

Co-authored-by: Denis Jelovina <[email protected]>
Co-authored-by: Vladimir Dimic <[email protected]>
# This is the commit message #189:

Resolve internal issue #561: Generalize ALP matrix over structure (#94)

* Generalize the implementation of ALP matrix over the structure.
* Conditionally enable square or rectangular constructors if the row and column dimensions are equal or not, respectively.
* Generalize the implementation of view_type within ALP matrix
* Highlight in a comment a minor issue with the compatibility check of matrix template arguments
# This is the commit message #190:

Resolve internal issue #563: Finish the implementation of select view (#97)

This commit adds missing pieces to enable select/permutation views over ALP matrix and ALP vector.
In addition, the interface now requires an ALP vector as the select/permutation vector.
# This is the commit message #191:

Resolve internal issue #572: Improve the design of factory for composed IMFs (#103)

This commit introduces the redesign of the factory for composed IMFs
The factory now exposes the type and offers a static create function that returns an object of the exposed type.
Removed separate composed_type type trait as its features are now implemented within the compose factory.

* Expand the unit test with the failing case that motivated this refactor.
* Fix the domain of the Zero IMF when creating a diagonal view.
# This is the commit message #192:

569 householder qr decomposition (#98)

QR decomposition real and complex version functional and smoketests enables

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #193:

559 dense symmetric tridiagonal eigensolver (#108)

-find zeroes of secular polynomial using bisections: function bisec_sec_eq
            this is not an optimal solution and it might be replaced later
-eigenvector calculations is not stable and gives non-ortohonal matrix for larger matrices
            an improved algorithm is required
-fold Matrix -> vector is missing and might be useful
-smoketest enables
-several commented issues should be addressed in the issue system 

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #194:

571 Generalize ALP Matrix and Vector over backend (#109)

* Restructure classes and structures related to ALP matrix
Add folder amf-based that serves as a base folder for AMF-based ALP containers.
* Minor code style fixes
* Remove obsolete documentation from base/matrix.hpp
* Refactor ALP vector and related stuff to be backend agnostic
* Specialize determine_poly_type to reference backend
* Reorder includes to be consistent with other files
* Revert exposing the type of underlying container
* Refactor functions to get band limits by separating logic to structure- and matrix-level
* Move functions querying matrix band limits to amf-base
* Move factory methods for constant matrices to amf-base
# This is the commit message #195:

missing review (#112)

type traits check

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #196:

Resolve internal issue #575: Fix the placement of operations and includes in blas files (#111)

* Move clear, resize, set, nnz and vector size-related getters to io.hpp
* Add missing io.hpp in alp root.
* Add base definitions in blas{0,1,2,3}.hpp

* Remove masked operation variants from blas1
* Remove obsolete operation definitions from blas{0,1,2,3} base
* Remove setElement taking a C++ scalar

* Fix the way NO_CAST_ASSERT and NO_CAST_OP_ASSERT are defined
* Add missing includes
* Fix code style issues
* Add asserts in base functions without an RC return value
# This is the commit message #197:

578 alp conjugate for vectors (#114)

- conjugate(vector) added to blas1.hpp
- unit test added

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #198:

574 trsv trsm via gauss jordan elimination (#116)

- trsv and trsm implemented trough backsubstitution  
- check numerical error
- conjugate  version for vectors
- use random_value() in checks
- add smoketest
- add matrix sizes checks

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #199:

571 Generalize internal containers over backend (#117)

* Move is_a to structures.hpp
* Add missing includes
* Generalize StorageBasedMatrix over backend
* Move friend functions of internal matrix-related classes to amf-base
* Move Functor-based matrix code to a separate file
* Move Storage-based matrix code to a separate file
* Remove unnecessary forward declarations
# This is the commit message #200:

571 Reorganize storage related classes (#119)

* Move storage-related types to amf-based folder
* Move storage-related type traits from matrix to storage files
* Move type traits that select storage implementation to storage file
# This is the commit message #201:

579 dense forward substitution (#118)

* linear system solver: forward substitution

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #202:

575 Move scalar version of set to IO (#120)

* Move scalar variant of set to IO
* Add missing undef for macro checking operator and parameter compatibility
# This is the commit message #203:

417 alp dense blocked cholesky integraton (#121)

- Cholesky versions:
        1) blocked  in-place
        2) blocked  out-of-place
        3) non-blocked in-place 
        4) non-blocked out-of-place
- Added structure checks
- forwardsubstitution inplace
- use Symmetric struct for reflector
- size checks
- test block size for corner cases
- added general view in the upper part of uppertriangular matrix

- Missing: SymmetricPositiveDefinite structure which is currently replace by Symmetric


Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #204:

582 remove getlength for alp vectors (#122)

* Replace calls to getLength with calls to size
* Move friend getLength to internal and refer to it from the public size function
# This is the commit message #205:

581 fix cmake infrastructure after rebasing on develop (#123)

* Clean up ALP-related cmake commands
* Clean up separation of include and include/graphblas cmake files
* Fix install path for transition
* Update the list of headers for alp backends
# This is the commit message #207:

symmetric eigensolver (#125)

- real-symmetric eigensolver
- depends on householder tridiagonalisation and tridiagonal eigensolver
                  Once done,  complex-hermitian version will use the same call

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #208:

586 dense cholesky for hermitian positive definite matrices (#126)

 Complex Cholesky - unblocked out-of-place
- enabled  cmake and smoketest
- add SymmetricPositiveDefinite and HermitianPositiveDefinite structures


Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #209:

560 alp reference sort primitive (#113)

* started adding sequential implementation of sort + unit test

* added vector iterators and internal begin/end

* begin/end logic defined privately

* Addressing first round of reviews

* Started definition of internal base classes for ALP relations.

* added first unit tests

* adapted sort to use the ALP relations' proposed draft

* Added missing basic relations and respective unit tests.

* setting (strict) poset requirement for sort's relation

* Addressing reviews.

* Addressing additional reviews: style fix, added grb-like set< use_index >, use alp::sort in symm_tridiag_dac_eigensolver algorithm.

* Code style fix

* Added review note in comment.

* Basing internal relation checks on internal ops apply (equal and not_equal only for now)
# This is the commit message #210:

585 dense symmetric hermitian positive definite matrix inversion (#127)

Positive definite matrix inversion
 - complex-hermitian
 - real-symmetric

Structures used SymmetricPositiveDefinite and HermitianPositiveDefinite

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #211:

590 dense lu factorisation (#128)

LU decomposition with pivoting (non-blocked version)

- Complex and real version
- Check all thin, square and flat input matrices
- Add Upper(Lower)Trapezoidal structures
- Enable smoketests

Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #212:

Refactoring includes

# This is the commit message #213:

596 add initial container and ops implementations for parallel dense backend (#129)

Introduce essential functionalities of the parallel backend to support allocating a matrix in parallel and setting its elements to a single value.

# This is the commit message #214:

601 complete support for replication across threads (#133)

* Replication factor in Distribution becomes static constexpr
* Take into account replication factor in thread grid layout
* Consider full thread grid coordinates (including rt) in containers and operations
* Add missing comments and fix a mistake in another comment
* Avoid calling distribution functions multiple times
* Pass thread coordinates using ThreadCoords object
* Compute number of threads within the distribution
* Encapsulate thread coordinates within the local coordinates
* Remove thread-grid related structures and getters since they are no longer needed
* Explain the reason behind hard-coded rt value
* Calculate block id inside distribution
* Rename distribution to a more descriptive name
# This is the commit message #215:

Refactoring includes (#137)


# This is the commit message #216:

592 dense svd factorisation (#132)

Singular Value Decomposition using Golub-Kahan algorithm

Bidiagonalisation of a general matrix (using Householder reflections) algorithm added.
Givens rotations added.
Checks for random matrices of size (n x n), (n x 2n) and (2n x n).
Implemented for real and complex matrices.
Numerical correctness is checked.
smoketests are enabled.
New structures are added:
 - LowerBidiagonal
 - UpperBidiagonal
 - RectangularDiagonal
and the following are updated:
- OrthogonalRows
- OrthogonalColumns 
TODO: 
 - Check Givens rotation and potential move it to the utils.
 - Adjust convergence threshold parameters and use relative convergence test.
 - Check for zeros in bidiagonal matrix and apply Givens rotation in order to eliminate an element in superdiagonal.


Co-authored-by: Denis Jelovina <[email protected]>
# This is the commit message #217:

Started drafting mxm using 2.5D algo

# This is the commit message #218:

WIP implementation of shifting + compute

# This is the commit message #219:

First complete draft

# This is the commit message #220:

compiling mxm test on general matrices

# This is the commit message #221:

shared mem (ge) mxm passing functional tests

# This is the commit message #222:

Tmp debugging info + fix in shift computation

# This is the commit message #223:

Refactoring includes

# This is the commit message #224:

- Fixing mismatching new/delete bug
- Enabling allocation/computation on sub-grid of threads
-  Added unit tests for non-cubic mxm

Author: hidanielesgit

CMakeLists.txt

@@ -59,6 +59,7 @@ if( WITH_HYPERDAGS_BACKEND )
 	endif()
 endif()
 option( WITH_NONBLOCKING_BACKEND "With Nonblocking backend" ON )
+option( WITH_ALP_OMP_BACKEND "With OMP Dense backend" ON )
 option( WITH_NUMA "With NUMA support" ON )
 option( LPF_INSTALL_PATH "Path to the LPF tools for the BSP1D and Hybrid backends" OFF )
 # the following options depend on LPF_INSTALL_PATH being set
@@ -307,6 +308,7 @@ set( WITH_REFERENCE_BACKEND_HEADERS OFF )
 set( WITH_ALP_REFERENCE_BACKEND_HEADERS OFF )
 set( WITH_OMP_BACKEND_HEADERS OFF )
 set( WITH_HYPERDAGS_BACKEND_HEADERS OFF )
+set( WITH_ALP_OMP_BACKEND_HEADERS OFF )
 
 # activate headers based on requested backends
 if( WITH_REFERENCE_BACKEND OR WITH_BSP1D_BACKEND OR WITH_NONBLOCKING_BACKEND )
@@ -327,6 +329,10 @@ if( WITH_ALP_REFERENCE_BACKEND )
 	set( WITH_ALP_REFERENCE_BACKEND_HEADERS ON )
 endif()
 
+if( WITH_ALP_OMP_BACKEND )
+	set( WITH_ALP_OMP_BACKEND_HEADERS ON )
+endif()
+
 add_subdirectory( include )
 
 ### BACKEND IMPLEMENTATIONS

bootstrap.sh

@@ -106,6 +106,7 @@ nonblocking=yes
 banshee=no
 lpf=no
 alp_reference=yes
+alp_omp=yes
 show=no
 FLAGS=$''
 LPF_INSTALL_PATH=
@@ -373,6 +374,11 @@ the current directory before invocation or confirm the deletion of its content w
 	else
 		CMAKE_OPTS+=" -DWITH_ALP_REFERENCE_BACKEND=ON"
 	fi
+	if [[ "${alp_omp}" == "no" ]]; then
+		CMAKE_OPTS+=" -DWITH_ALP_OMP_BACKEND=OFF"
+	else
+		CMAKE_OPTS+=" -DWITH_ALP_OMP_BACKEND=ON"
+	fi
 	if [[ "${lpf}" == "yes" ]]; then
 		CMAKE_OPTS+=" -DLPF_INSTALL_PATH='${ABSOLUTE_LPF_INSTALL_PATH}'"
 	fi

cmake/AddGRBInstall.cmake

@@ -46,6 +46,7 @@ set( ALP_UTILS_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}" )
 set( SHMEM_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/sequential" )
 set( HYPERDAGS_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/hyperdags" )
 set( ALP_REFERENCE_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/alp/reference" )
+set( ALP_OMP_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/alp/omp" )
 set( BSP1D_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/spmd" )
 set( HYBRID_BACKEND_INSTALL_DIR "${BINARY_LIBRARIES_INSTALL_DIR}/hybrid" )
 
@@ -154,6 +155,14 @@ if( WITH_NONBLOCKING_BACKEND )
 	)
 endif()
 
+if( WITH_ALP_OMP_BACKEND )
+	addBackendWrapperGenOptions( "alp_omp"
+		COMPILE_DEFINITIONS "${ALP_OMP_SELECTION_DEFS}"
+		LINK_FLAGS "'${ALP_OMP_BACKEND_INSTALL_DIR}/lib${BACKEND_LIBRARY_OUTPUT_NAME}.a'"
+			"'${ALP_UTILS_INSTALL_DIR}/lib${ALP_UTILS_LIBRARY_OUTPUT_NAME}.a'" "${NUMA_LFLAG}"
+	)
+endif()
+
 # distributed memory backends
 if( WITH_BSP1D_BACKEND OR WITH_HYBRID_BACKEND )
 	assert_valid_variables( LPFRUN LPFCPP )

cmake/AddGRBVars.cmake

@@ -34,6 +34,7 @@ set( HYBRID_BACKEND_DEFAULT_NAME "backend_hybrid" )
 set( HYPERDAGS_BACKEND_DEFAULT_NAME "backend_hyperdags" )
 set( NONBLOCKING_BACKEND_DEFAULT_NAME "backend_nonblocking" )
 set( ALP_REFERENCE_BACKEND_DEFAULT_NAME "backend_alp_reference" )
+set( ALP_OMP_BACKEND_DEFAULT_NAME "backend_alp_omp" )
 
 
 ### COMPILER DEFINITIONS FOR HEADERS INCLUSION AND FOR BACKEND SELECTION
@@ -45,6 +46,7 @@ set( HYPERDAGS_INCLUDE_DEFS "_GRB_WITH_HYPERDAGS" )
 set( NONBLOCKING_INCLUDE_DEFS "_GRB_WITH_NONBLOCKING" )
 set( LPF_INCLUDE_DEFS "_GRB_WITH_LPF" )
 set( ALP_REFERENCE_INCLUDE_DEFS "_ALP_WITH_REFERENCE" )
+set( ALP_OMP_INCLUDE_DEFS "_ALP_WITH_OMP;_ALP_OMP_WITH_REFERENCE" )
 
 # compiler definitions to select a backend
 set( REFERENCE_SELECTION_DEFS "_GRB_BACKEND=reference" )
@@ -55,6 +57,10 @@ set( HYPERDAGS_SELECTION_DEFS
 )
 set( NONBLOCKING_SELECTION_DEFS "_GRB_BACKEND=nonblocking" )
 set( ALP_REFERENCE_SELECTION_DEFS "_ALP_BACKEND=reference" )
+set( ALP_OMP_SELECTION_DEFS
+		"_ALP_BACKEND=omp"
+		"_ALP_SECONDARY_BACKEND=reference"
+)
 set( BSP1D_SELECTION_DEFS
 		"_GRB_BACKEND=BSP1D"
 		"_GRB_BSP1D_BACKEND=reference"
@@ -68,7 +74,7 @@ set( HYBRID_SELECTION_DEFS
 set( NO_NUMA_DEF "_GRB_NO_LIBNUMA" )
 
 ### **ALL** BACKENDS, EVEN IF NOT ENABLED BY USER
-set( ALL_BACKENDS "reference" "reference_omp" "hyperdags" "nonblocking" "bsp1d" "hybrid" "alp_reference" )
+set( ALL_BACKENDS "reference" "reference_omp" "hyperdags" "nonblocking" "bsp1d" "hybrid" "alp_reference" "alp_omp" )
 
 # list of user-enabled backends, for tests and wrapper scripts (do not change!)
 set( AVAILABLE_BACKENDS "" )
@@ -98,6 +104,10 @@ if( WITH_ALP_REFERENCE_BACKEND )
 	list( APPEND AVAILABLE_BACKENDS "alp_reference" )
 endif()
 
+if( WITH_ALP_OMP_BACKEND )
+	list( APPEND AVAILABLE_BACKENDS "alp_omp" )
+endif()
+
 # distributed memory backends
 if( WITH_BSP1D_BACKEND )
 	list( APPEND AVAILABLE_BACKENDS "bsp1d" )

include/CMakeLists.txt

@@ -108,7 +108,9 @@ install( TARGETS alp_utils_headers EXPORT GraphBLASTargets
 	INCLUDES DESTINATION "${INCLUDE_INSTALL_DIR}"
 )
 
-if( WITH_ALP_REFERENCE_BACKEND_HEADERS )
+if( WITH_ALP_REFERENCE_BACKEND_HEADERS OR
+	WITH_ALP_OMP_BACKEND_HEADERS
+)
 	add_subdirectory( alp )
 endif()
 

include/alp.hpp

@@ -134,7 +134,7 @@
 // clang-format re-ordering
 #if 1
 // load active configuration
-// #include <alp/config.hpp> //defines _ALP_BACKEND and _WITH_BSP
+ #include <alp/config.hpp> //defines _ALP_BACKEND
 #endif
 
 // #pragma message "Included ALP.hpp"

include/alp/CMakeLists.txt

@@ -22,6 +22,7 @@
 # set a default backend (if they want to do so).
 #
 assert_defined_variables( ALP_REFERENCE_INCLUDE_DEFS WITH_ALP_REFERENCE_BACKEND_HEADERS )
+assert_defined_variables( ALP_OMP_INCLUDE_DEFS WITH_ALP_OMP_BACKEND_HEADERS )
 assert_valid_variables( INCLUDE_INSTALL_DIR )
 
 # to avoid flaky acrobatics with regex or glob expressions, copy main files directly
@@ -91,6 +92,22 @@ if( WITH_ALP_REFERENCE_BACKEND )
 	)
 endif()
 
+if( WITH_ALP_OMP_BACKEND_HEADERS )
+	add_library( backend_alp_omp_headers INTERFACE )
+	target_link_libraries( backend_alp_omp_headers INTERFACE backend_headers_nodefs )
+	target_link_libraries( backend_alp_omp_headers INTERFACE OpenMP::OpenMP_CXX )
+	target_compile_definitions( backend_alp_omp_headers INTERFACE "${ALP_OMP_INCLUDE_DEFS}" )
+
+	install( TARGETS backend_alp_omp_headers EXPORT GraphBLASTargets )
+endif()
+
+if( WITH_ALP_OMP_BACKEND )
+	install( DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}/alp/omp"
+		DESTINATION "${ALP_INCLUDE_INSTALL_DIR}/alp_omp"
+		FILES_MATCHING REGEX "${HEADERS_REGEX}"
+	)
+endif()
+
 # this target lists the algorithms implemented on top of the generic functionalities,
 # hence it depends only on backend_headers_nodefs
 add_library( alp_algorithms INTERFACE )

include/alp/algorithms/householder_bidiag.hpp

@@ -0,0 +1,204 @@
+/*
+ *   Copyright 2021 Huawei Technologies Co., Ltd.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include <iostream>
+#include <sstream>
+
+#include <alp.hpp>
+#include <graphblas/utils/iscomplex.hpp> // use from grb
+#ifdef DEBUG
+#include "../tests/utils/print_alp_containers.hpp"
+#endif
+
+namespace alp {
+
+	namespace algorithms {
+
+		/**
+		 * Given a general matrix H perform an inplace Householder reflections in order to
+		 * eliminate column elements H[i+d:,i] (below diagonal or subdigonal), H = U H.
+		 *
+		 * @tparam D        Data element type
+		 * @tparam Ring     Type of the semiring used in the computation
+		 * @tparam Minus    Type minus operator used in the computation
+		 * @tparam Divide   Type of divide operator used in the computation
+		 * @param[in,out]   U updated orthogonal matrix
+		 * @param[in,out]   H updated general matrix with column
+		 *                         elements H[i+d:,i] eliminated
+		 * @param[in]       i column to eliminate
+		 * @param[in]       d offset from diagonal to eliminate, default 0
+		 * @param[in]       ring A semiring for operations
+		 * @return RC       SUCCESS if the execution was correct
+		 */
+		template<
+			typename MatH,
+			typename MatU,
+			typename IndexType,
+			typename D = typename MatH::value_type,
+			class Ring = Semiring< operators::add< D >, operators::mul< D >, identities::zero, identities::one >,
+			class Minus = operators::subtract< D >,
+			class Divide = operators::divide< D >,
+			std::enable_if_t<
+				structures::is_a< typename MatH::structure, structures::General >::value &&
+				structures::is_a< typename MatU::structure, structures::Orthogonal >::value &&
+				is_semiring< Ring >::value &&
+				is_operator< Minus >::value &&
+				is_operator< Divide >::value
+			> * = nullptr
+		>
+		RC elminate_below_ith_diag(
+			const IndexType i,
+			MatH &H,
+			MatU &U,
+			const IndexType d = 0,
+			const Ring &ring = Ring(),
+			const Minus &minus = Minus(),
+			const Divide &divide = Divide()
+		) {
+			RC rc = SUCCESS;
+
+			const Scalar< D > zero( ring.template getZero< D >() );
+			const Scalar< D > one( ring.template getOne< D >() );
+
+			const IndexType m = nrows( H );
+			const IndexType n = ncols( H );
+
+			// v=copy(A0[i+d:,i])
+			auto a = get_view( H, utils::range( i + d, m ), i );
+			Vector< D > v( m - ( i + d ) );
+			rc = rc ? rc : alp::set( v, a );
+
+			// alpha=v[0]/abs(v[0])
+			Scalar< D > alpha( zero );
+			auto v0 = get_view( v, utils::range( 0, 1 ) );
+			rc = rc ? rc : foldl( alpha, v0, ring.getAdditiveMonoid() );
+			rc = rc ? rc : foldl( alpha, Scalar< D >( std::abs( *alpha ) ), divide );
+
+			// alpha=alpha*norm(v)
+			Scalar< D > norm_v1( zero );
+			rc = rc ? rc : norm2( norm_v1, v, ring );
+			rc = rc ? rc : foldl( alpha, norm_v1, ring.getMultiplicativeOperator() );
+
+			// v[0]=v[0]-alpha
+			rc = rc ? rc : foldl( v0, alpha, minus );
+
+			// v=v/norm(v)
+			Scalar< D > norm_v2( zero );
+			rc = rc ? rc : norm2( norm_v2, v, ring );
+			rc = rc ? rc : foldl( v, norm_v2, divide );
+
+			//P1=zeros((m-(i+d),m-(i+d))).astype(complex)
+			//P1=P1-2*outer(v,conjugate(v))
+			auto vvh = outer( v, ring.getMultiplicativeOperator() );
+			Matrix< D, typename decltype( vvh )::structure, Dense > reflector( m - ( i + d ) );
+			rc = rc ? rc : alp::set( reflector, vvh );
+			rc = rc ? rc : foldl( reflector, Scalar< D > ( -2 ), ring.getMultiplicativeOperator() );
+
+			// A0=P.dot(A0)
+			auto Hupdate = get_view( H, utils::range( i + d, m ), utils::range( 0, n ) );
+			Matrix< D, structures::General, Dense > Temp1( m - ( i + d ) , n );
+			rc = rc ? rc : alp::set( Temp1, Hupdate );
+			rc = rc ? rc : mxm( Hupdate, reflector, Temp1, ring );
+
+			// Uk=Uk.dot(P)
+			auto Uupdate = get_view< structures::OrthogonalColumns >( U, utils::range( 0, m ), utils::range( i + d, m ) );
+			Matrix< D, structures::OrthogonalColumns, Dense > Temp2( m, m - ( i + d ) );
+			rc = rc ? rc : alp::set( Temp2, Uupdate );
+			rc = rc ? rc : mxm( Uupdate, Temp2, reflector, ring );
+
+			return rc;
+		}
+
+		/**
+		 * Computes Householder (inplace) bidiagonalisation of general matrix \f$H = U B V \f$
+		 * where \a H is general (complex or real),
+		 * \a U orthogonal, \a B is bidiagonal and  \a V orthogonal.
+		 *
+		 * @tparam D        Data element type
+		 * @tparam Ring     Type of the semiring used in the computation
+		 * @tparam Minus    Type minus operator used in the computation
+		 * @tparam Divide   Type of divide operator used in the computation
+		 * @param[in,out]   U updated orthogonal matrix
+		 * @param[in,out]   V updated orthogonal matrix
+		 * @param[in,out]   H input general matrix, output bidiagonal matrix (B)
+		 * @param[in]       ring A semiring for operations
+		 * @return RC       SUCCESS if the execution was correct
+		 *
+		 */
+		template<
+			typename MatH,
+			typename D = typename MatH::value_type,
+			typename MatU,
+			typename MatV,
+			class Ring = Semiring< operators::add< D >, operators::mul< D >, identities::zero, identities::one >,
+			class Minus = operators::subtract< D >,
+			class Divide = operators::divide< D >,
+			std::enable_if_t<
+				is_matrix< MatH >::value &&
+				is_matrix< MatU >::value &&
+				is_matrix< MatV >::value &&
+				structures::is_a< typename MatH::structure, structures::General >::value &&
+				structures::is_a< typename MatU::structure, structures::Orthogonal >::value &&
+				structures::is_a< typename MatV::structure, structures::Orthogonal >::value &&
+				is_semiring< Ring >::value &&
+				is_operator< Minus >::value &&
+				is_operator< Divide >::value
+			> * = nullptr
+		>
+		RC householder_bidiag(
+			MatU &U,
+			MatH &H,
+			MatV &V,
+			const Ring &ring = Ring(),
+			const Minus &minus = Minus(),
+			const Divide &divide = Divide()
+		) {
+			RC rc = SUCCESS;
+
+			const Scalar< D > zero( ring.template getZero< D >() );
+			const Scalar< D > one( ring.template getOne< D >() );
+
+			const size_t m = nrows( H );
+			const size_t n = ncols( H );
+
+			// check sizes
+			if(
+				( ncols( U ) != nrows( H ) ) ||
+				( ncols( H ) != nrows( V ) )
+			) {
+				std::cerr << "Incompatible sizes in householder_bidiag.\n";
+				return FAILED;
+			}
+
+			//for i in range(min(n,m)):
+			for( size_t i = 0; i < std::min( n, m ); ++i ) {
+				// eliminate column elements below ith diagonal element
+				if( i < std::min( n, m - 1 ) ) {
+					rc = rc ? rc : elminate_below_ith_diag( i, H, U, static_cast< size_t >( 0 ), ring, minus, divide );
+				}
+				// eliminate row elements to the right from (i+1)th diagonal element
+				if( i < std::min( n - 2, m ) ) {
+					auto HT = get_view< alp::view::transpose >( H );
+					auto VT = get_view< alp::view::transpose >( V );
+					rc = rc ? rc : elminate_below_ith_diag( i, HT, VT, static_cast< size_t >( 1 ), ring, minus, divide );
+				}
+			}
+
+			return rc;
+
+		}
+	} // namespace algorithms
+} // namespace alp

include/alp/algorithms/householder_lu.hpp

@@ -96,8 +96,10 @@ namespace alp {
 				( nrows( U ) != kk ) ||
 				( ncols( L ) != kk )
 			) {
+#ifdef DEBUG
 				std::cerr << " n, kk, m = " << n << ", "  << kk << ", " << m << "\n";
 				std::cerr << "Incompatible sizes in householder_lu.\n";
+#endif
 				return FAILED;
 			}