Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly
Sign up
Appearance settings

Commit 3369013

Browse files
authored
Sparse algebra support with OOP API (#760)
2 parents 427bc68 + 22cd23e commit 3369013

14 files changed

+3079
-0
lines changed

‎doc/specs/stdlib_sparse.md‎

Lines changed: 364 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,364 @@
1+
---
2+
title: sparse
3+
---
4+
5+
# The `stdlib_sparse` module
6+
7+
[TOC]
8+
9+
## Introduction
10+
11+
The `stdlib_sparse` module provides derived types for standard sparse matrix data structures. It also provides math kernels such as sparse matrix-vector product and conversion between matrix types.
12+
13+
## Sparse matrix derived types
14+
15+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
16+
### The `sparse_type` abstract derived type
17+
#### Status
18+
19+
Experimental
20+
21+
#### Description
22+
The parent `sparse_type` is as an abstract derived type holding the basic common meta data needed to define a sparse matrix, as well as shared APIs. All sparse matrix flavors are extended from the `sparse_type`.
23+
24+
```Fortran
25+
type, public, abstract :: sparse_type
26+
integer :: nrows !! number of rows
27+
integer :: ncols !! number of columns
28+
integer :: nnz !! number of non-zero values
29+
integer :: storage !! assumed storage symmetry
30+
end type
31+
```
32+
33+
The storage integer label should be assigned from the module's internal enumerator containing the following three enums:
34+
35+
```Fortran
36+
enum, bind(C)
37+
enumerator :: sparse_full !! Full Sparse matrix (no symmetry considerations)
38+
enumerator :: sparse_lower !! Symmetric Sparse matrix with triangular inferior storage
39+
enumerator :: sparse_upper !! Symmetric Sparse matrix with triangular supperior storage
40+
end enum
41+
```
42+
In the following, all sparse kinds will be presented in two main flavors: a data-less type `<matrix>_type` useful for topological graph operations. And real/complex valued types `<matrix>_<kind>_type` containing the `data` buffer for the matrix values. The following rectangular matrix will be used to showcase how each sparse matrix holds the data internally:
43+
44+
$$ M = \begin{bmatrix}
45+
9 & 0 & 0 & 0 & -3 \\
46+
4 & 7 & 0 & 0 & 0 \\
47+
0 & 8 & -1 & 8 & 0 \\
48+
4 & 0 & 5 & 6 & 0 \\
49+
\end{bmatrix} $$
50+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
51+
### `COO`: The COOrdinates compressed sparse format
52+
#### Status
53+
54+
Experimental
55+
56+
#### Description
57+
The `COO`, triplet or `ijv` format defines all non-zero elements of the matrix by explicitly allocating the `i,j` index and the value of the matrix. While some implementations use separate `row` and `col` arrays for the index, here we use a 2D array in order to promote fast memory acces to `ij`.
58+
59+
```Fortran
60+
type(COO_sp_type) :: COO
61+
call COO%malloc(4,5,10)
62+
COO%data(:) = real([9,-3,4,7,8,-1,8,4,5,6])
63+
COO%index(1:2,1) = [1,1]
64+
COO%index(1:2,2) = [1,5]
65+
COO%index(1:2,3) = [2,1]
66+
COO%index(1:2,4) = [2,2]
67+
COO%index(1:2,5) = [3,2]
68+
COO%index(1:2,6) = [3,3]
69+
COO%index(1:2,7) = [3,4]
70+
COO%index(1:2,8) = [4,1]
71+
COO%index(1:2,9) = [4,3]
72+
COO%index(1:2,10) = [4,4]
73+
```
74+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
75+
### `CSR`: The Compressed Sparse Row or Yale format
76+
#### Status
77+
78+
Experimental
79+
80+
#### Description
81+
The Compressed Sparse Row or Yale format `CSR` stores the matrix structure by compressing the row indices with a counter pointer `rowptr` enabling to know the first and last non-zero column index `col` of the given row.
82+
83+
```Fortran
84+
type(CSR_sp_type) :: CSR
85+
call CSR%malloc(4,5,10)
86+
CSR%data(:) = real([9,-3,4,7,8,-1,8,4,5,6])
87+
CSR%col(:) = [1,5,1,2,2,3,4,1,3,4]
88+
CSR%rowptr(:) = [1,3,5,8,11]
89+
```
90+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
91+
### `CSC`: The Compressed Sparse Column format
92+
#### Status
93+
94+
Experimental
95+
96+
#### Description
97+
The Compressed Sparse Colum `CSC` is similar to the `CSR` format but values are accesed first by column, thus an index counter is given by `colptr` which enables to know the first and last non-zero row index of a given colum.
98+
99+
```Fortran
100+
type(CSC_sp_type) :: CSC
101+
call CSC%malloc(4,5,10)
102+
CSC%data(:) = real([9,4,4,7,8,-1,5,8,6,-3])
103+
CSC%row(:) = [1,2,4,2,3,3,4,3,4,1]
104+
CSC%colptr(:) = [1,4,6,8,10,11]
105+
```
106+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
107+
### `ELLPACK`: ELL-pack storage format
108+
#### Status
109+
110+
Experimental
111+
112+
#### Description
113+
The `ELL` format stores data in a dense matrix of $nrows \times K$ in column major order. By imposing a constant number of elements per row $K,ドル this format will incur in additional zeros being stored, but it enables efficient vectorization as memory acces is carried out by constant sized strides.
114+
115+
```Fortran
116+
type(ELL_sp_type) :: ELL
117+
call ELL%malloc(num_rows=4,num_cols=5,num_nz_row=3)
118+
ELL%data(1,1:3) = real([9,-3,0])
119+
ELL%data(2,1:3) = real([4,7,0])
120+
ELL%data(3,1:3) = real([8,-1,8])
121+
ELL%data(4,1:3) = real([4,5,6])
122+
123+
ELL%index(1,1:3) = [1,5,0]
124+
ELL%index(2,1:3) = [1,2,0]
125+
ELL%index(3,1:3) = [2,3,4]
126+
ELL%index(4,1:3) = [1,3,4]
127+
```
128+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
129+
### `SELL-C`: The Sliced ELLPACK with Constant blocks format
130+
#### Status
131+
132+
Experimental
133+
134+
#### Description
135+
The Sliced ELLPACK format `SELLC` is a variation of the `ELLPACK` format. This modification reduces the storage size compared to the `ELLPACK` format but maintaining its efficient data access scheme. It can be seen as an intermediate format between `CSR` and `ELLPACK`. For more details read [the reference](https://arxiv.org/pdf/1307.6209v1)
136+
137+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
138+
## `add`- sparse matrix data accessors
139+
140+
### Status
141+
142+
Experimental
143+
144+
### Description
145+
Type-bound procedures to enable adding data in a sparse matrix.
146+
147+
### Syntax
148+
149+
`call matrix%add(i,j,v)` or
150+
`call matrix%add(i(:),j(:),v(:,:))`
151+
152+
### Arguments
153+
154+
`i`: Shall be an integer value or rank-1 array. It is an `intent(in)` argument.
155+
156+
`j`: Shall be an integer value or rank-1 array. It is an `intent(in)` argument.
157+
158+
`v`: Shall be a `real` or `complex` value or rank-2 array. The type shall be in accordance to the declared sparse matrix object. It is an `intent(in)` argument.
159+
160+
## `at`- sparse matrix data accessors
161+
162+
### Status
163+
164+
Experimental
165+
166+
### Description
167+
Type-bound procedures to enable requesting data from a sparse matrix.
168+
169+
### Syntax
170+
171+
`v = matrix%at(i,j)`
172+
173+
### Arguments
174+
175+
`i` : Shall be an integer value. It is an `intent(in)` argument.
176+
177+
`j` : Shall be an integer value. It is an `intent(in)` argument.
178+
179+
`v` : Shall be a `real` or `complex` value in accordance to the declared sparse matrix object. If the `ij` tuple is within the sparse pattern, `v` contains the value in the data buffer. If the `ij` tuple is outside the sparse pattern, `v` is equal `0`. If the `ij` tuple is outside the matrix pattern `(nrows,ncols)`, `v` is `NaN`.
180+
181+
## Example
182+
```fortran
183+
{!example/linalg/example_sparse_data_accessors.f90!}
184+
```
185+
186+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
187+
## `spmv` - Sparse Matrix-Vector product
188+
189+
### Status
190+
191+
Experimental
192+
193+
### Description
194+
195+
Provide sparse matrix-vector product kernels for the current supported sparse matrix types.
196+
197+
$$y=\alpha*op(M)*x+\beta*y$$
198+
199+
### Syntax
200+
201+
`call ` [[stdlib_sparse_spmv(module):spmv(interface)]] `(matrix,vec_x,vec_y [,alpha,beta,op])`
202+
203+
### Arguments
204+
205+
`matrix`: Shall be a `real` or `complex` sparse type matrix. It is an `intent(in)` argument.
206+
207+
`vec_x`: Shall be a rank-1 or rank-2 array of `real` or `complex` type array. It is an `intent(in)` argument.
208+
209+
`vec_y`: Shall be a rank-1 or rank-2 array of `real` or `complex` type array. . It is an `intent(inout)` argument.
210+
211+
`alpha`, `optional` : Shall be a scalar value of the same type as `vec_x`. Default value `alpha=1`. It is an `intent(in)` argument.
212+
213+
`beta`, `optional` : Shall be a scalar value of the same type as `vec_x`. Default value `beta=0`. It is an `intent(in)` argument.
214+
215+
`op`, `optional`: In-place operator identifier. Shall be a `character(1)` argument. It can have any of the following values: `N`: no transpose, `T`: transpose, `H`: hermitian or complex transpose. These values are provided as constants by the `stdlib_sparse` module: `sparse_op_none`, `sparse_op_transpose`, `sparse_op_hermitian`
216+
217+
<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
218+
## Sparse matrix to matrix conversions
219+
220+
### Status
221+
222+
Experimental
223+
224+
### Description
225+
226+
This module provides facility functions for converting between storage formats.
227+
228+
### Syntax
229+
230+
`call ` [[stdlib_sparse_conversion(module):coo2ordered(interface)]] `(coo[,sort_data])`
231+
232+
### Arguments
233+
234+
`COO` : Shall be any `COO` type. The same object will be returned with the arrays reallocated to the correct size after removing duplicates. It is an `intent(inout)` argument.
235+
236+
`sort_data`, `optional` : Shall be a `logical` argument to determine whether data in the COO graph should be sorted while sorting the index array, default `.false.`. It is an `intent(in)` argument.
237+
238+
### Syntax
239+
240+
`call ` [[stdlib_sparse_conversion(module):from_ijv(interface)]] `(sparse,row,col[,data,nrows,ncols,num_nz_rows,chunk])`
241+
242+
### Arguments
243+
244+
`sparse` : Shall be a `COO`, `CSR`, `ELL` or `SELLC` type. The graph object will be returned with a canonical shape after sorting and removing duplicates from the `(row,col,data)` triplet. If the graph is `COO_type` no data buffer is allowed. It is an `intent(inout)` argument.
245+
246+
`row` : rows index array. It is an `intent(in)` argument.
247+
248+
`col` : columns index array. It is an `intent(in)` argument.
249+
250+
`data`, `optional`: `real` or `complex` data array. It is an `intent(in)` argument.
251+
252+
`nrows`, `optional`: number of rows, if not given it will be computed from the `row` array. It is an `intent(in)` argument.
253+
254+
`ncols`, `optional`: number of columns, if not given it will be computed from the `col` array. It is an `intent(in)` argument.
255+
256+
`num_nz_rows`, `optional`: number of non zeros per row, only valid in the case of an `ELL` matrix, by default it will computed from the largest row. It is an `intent(in)` argument.
257+
258+
`chunk`, `optional`: chunk size, only valid in the case of a `SELLC` matrix, by default it will be taken from the `SELLC` default attribute chunk size. It is an `intent(in)` argument.
259+
260+
## Example
261+
```fortran
262+
{!example/linalg/example_sparse_from_ijv.f90!}
263+
```
264+
### Syntax
265+
266+
`call ` [[stdlib_sparse_conversion(module):diag(interface)]] `(matrix,diagonal)`
267+
268+
### Arguments
269+
270+
`matrix` : Shall be a `dense`, `COO`, `CSR` or `ELL` type. It is an `intent(in)` argument.
271+
272+
`diagonal` : A rank-1 array of the same type as the `matrix`. It is an `intent(inout)` and `allocatable` argument.
273+
274+
#### Note
275+
If the `diagonal` array has not been previously allocated, the `diag` subroutine will allocate it using the `nrows` of the `matrix`.
276+
277+
### Syntax
278+
279+
`call ` [[stdlib_sparse_conversion(module):dense2coo(interface)]] `(dense,coo)`
280+
281+
### Arguments
282+
283+
`dense` : Shall be a rank-2 array of `real` or `complex` type. It is an `intent(in)` argument.
284+
285+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(out)` argument.
286+
287+
### Syntax
288+
289+
`call ` [[stdlib_sparse_conversion(module):coo2dense(interface)]] `(coo,dense)`
290+
291+
### Arguments
292+
293+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(in)` argument.
294+
295+
`dense` : Shall be a rank-2 array of `real` or `complex` type. It is an `intent(out)` argument.
296+
297+
### Syntax
298+
299+
`call ` [[stdlib_sparse_conversion(module):coo2csr(interface)]] `(coo,csr)`
300+
301+
### Arguments
302+
303+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(in)` argument.
304+
305+
`csr` : Shall be a `CSR` type of `real` or `complex` type. It is an `intent(out)` argument.
306+
307+
### Syntax
308+
309+
`call ` [[stdlib_sparse_conversion(module):coo2csc(interface)]] `(coo,csc)`
310+
311+
### Arguments
312+
313+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(in)` argument.
314+
315+
`csc` : Shall be a `CSC` type of `real` or `complex` type. It is an `intent(out)` argument.
316+
317+
### Syntax
318+
319+
`call ` [[stdlib_sparse_conversion(module):csr2coo(interface)]] `(csr,coo)`
320+
321+
### Arguments
322+
323+
`csr` : Shall be a `CSR` type of `real` or `complex` type. It is an `intent(in)` argument.
324+
325+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(out)` argument.
326+
327+
### Syntax
328+
329+
`call ` [[stdlib_sparse_conversion(module):csr2sellc(interface)]] `(csr,sellc[,chunk])`
330+
331+
### Arguments
332+
333+
`csr` : Shall be a `CSR` type of `real` or `complex` type. It is an `intent(in)` argument.
334+
335+
`sellc` : Shall be a `SELLC` type of `real` or `complex` type. It is an `intent(out)` argument.
336+
337+
`chunk`, `optional`: chunk size for the Sliced ELLPACK format. It is an `intent(in)` argument.
338+
339+
### Syntax
340+
341+
`call ` [[stdlib_sparse_conversion(module):csr2sellc(interface)]] `(csr,ell[,num_nz_rows])`
342+
343+
### Arguments
344+
345+
`csr` : Shall be a `CSR` type of `real` or `complex` type. It is an `intent(in)` argument.
346+
347+
`ell` : Shall be a `ELL` type of `real` or `complex` type. It is an `intent(out)` argument.
348+
349+
`num_nz_rows`, `optional`: number of non zeros per row. If not give, it will correspond to the size of the longest row in the `CSR` matrix. It is an `intent(in)` argument.
350+
351+
### Syntax
352+
353+
`call ` [[stdlib_sparse_conversion(module):csc2coo(interface)]] `(csc,coo)`
354+
355+
### Arguments
356+
357+
`csc` : Shall be a `CSC` type of `real` or `complex` type. It is an `intent(in)` argument.
358+
359+
`coo` : Shall be a `COO` type of `real` or `complex` type. It is an `intent(out)` argument.
360+
361+
## Example
362+
```fortran
363+
{!example/linalg/example_sparse_spmv.f90!}
364+
```

‎example/linalg/CMakeLists.txt‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,9 @@ ADD_EXAMPLE(get_norm)
3333
ADD_EXAMPLE(solve1)
3434
ADD_EXAMPLE(solve2)
3535
ADD_EXAMPLE(solve3)
36+
ADD_EXAMPLE(sparse_from_ijv)
37+
ADD_EXAMPLE(sparse_data_accessors)
38+
ADD_EXAMPLE(sparse_spmv)
3639
ADD_EXAMPLE(svd)
3740
ADD_EXAMPLE(svdvals)
3841
ADD_EXAMPLE(determinant)

0 commit comments

Comments
(0)

AltStyle によって変換されたページ (->オリジナル) /