General sparsity class. More...
#include <sparsity.hpp>
The storage format is a compressed column storage (CCS) format.
In this format, the structural non-zero elements are stored in column-major order, starting from the upper left corner of the matrix and ending in the lower right corner.
In addition to the dimension (size1(), size2()), (i.e. the number of rows and the number of columns respectively), there are also two vectors of integers:
Note that with this format, it is cheap to loop over all the non-zero elements of a particular column, at constant time per element, but expensive to jump to access a location (i, j).
If the matrix is dense, i.e. length(row) == size1()*size2(), the format reduces to standard dense column major format, which allows access to an arbitrary element in constant time.
Since the object is reference counted (it inherits from SharedObject), several matrices are allowed to share the same sparsity pattern.
The implementations of methods marked as such in this class has been taken from the CSparse package and modified to fit CasADi data structures and separation of sparsity pattern calculation and numerical evaluation. These functions are Copyright(c) Timothy A. Davis, 2006-2009 and licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_b9
Definition at line 96 of file sparsity.hpp.
Public Member Functions | |
Sparsity (casadi_int dummy=0) | |
Default constructor. More... | |
Sparsity (casadi_int nrow, casadi_int ncol) | |
Pattern with all structural zeros. More... | |
Sparsity (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &colind, const std::vector< casadi_int > &row, bool order_rows=false) | |
Construct from sparsity pattern vectors given in compressed column storage format. More... | |
Sparsity (const std::pair< casadi_int, casadi_int > &rc) | |
Create a sparse matrix with all structural zeros. More... | |
const std::vector< casadi_int > | permutation_vector (bool invert=false) const |
Construct permutation vector from permutation matrix. More... | |
Sparsity | get_diag (std::vector< casadi_int > &mapping) const |
std::vector< casadi_int > | compress () const |
Compress a sparsity pattern. More... | |
bool | operator!= (const Sparsity &y) const |
Check if two sparsity patterns are difference. More... | |
bool | is_stacked (const Sparsity &y, casadi_int n) const |
Check if pattern is horizontal repeat of another. More... | |
Dict | info () const |
void | to_file (const std::string &filename, const std::string &format_hint="") const |
std::string | serialize () const |
Serialize. More... | |
void | serialize (SerializingStream &s) const |
Serialize an object. More... | |
std::vector< casadi_int > | get_row () const |
Get the row for each non-zero entry. More... | |
std::vector< casadi_int > | get_colind () const |
Get the column index for each column. More... | |
casadi_int | colind (casadi_int cc) const |
Get a reference to the colindex of column cc (see class description) More... | |
casadi_int | row (casadi_int el) const |
Get the row of a non-zero element. More... | |
std::vector< casadi_int > | get_col () const |
Get the column for each non-zero entry. More... | |
void | resize (casadi_int nrow, casadi_int ncol) |
Resize. More... | |
casadi_int | add_nz (casadi_int rr, casadi_int cc) |
Get the index of a non-zero element. More... | |
casadi_int | get_nz (casadi_int rr, casadi_int cc) const |
Get the index of an existing non-zero element. More... | |
bool | has_nz (casadi_int rr, casadi_int cc) const |
Returns true if the pattern has a non-zero at location rr, cc. More... | |
std::vector< casadi_int > | get_nz (const std::vector< casadi_int > &rr, const std::vector< casadi_int > &cc) const |
Get a set of non-zero element. More... | |
void | get_nz (std::vector< casadi_int > &indices) const |
Get the nonzero index for a set of elements. More... | |
std::vector< casadi_int > | get_lower () const |
Get nonzeros in lower triangular part. More... | |
std::vector< casadi_int > | get_upper () const |
Get nonzeros in upper triangular part. More... | |
void | get_ccs (std::vector< casadi_int > &colind, std::vector< casadi_int > &row) const |
Get the sparsity in compressed column storage (CCS) format. More... | |
void | get_crs (std::vector< casadi_int > &rowind, std::vector< casadi_int > &col) const |
Get the sparsity in compressed row storage (CRS) format. More... | |
void | get_triplet (std::vector< casadi_int > &row, std::vector< casadi_int > &col) const |
Get the sparsity in sparse triplet format. More... | |
Sparsity | sub (const std::vector< casadi_int > &rr, const std::vector< casadi_int > &cc, std::vector< casadi_int > &mapping, bool ind1=false) const |
Get a submatrix. More... | |
Sparsity | sub (const std::vector< casadi_int > &rr, const Sparsity &sp, std::vector< casadi_int > &mapping, bool ind1=false) const |
Get a set of elements. More... | |
Sparsity | T () const |
Transpose the matrix. More... | |
Sparsity | transpose (std::vector< casadi_int > &mapping, bool invert_mapping=false) const |
Transpose the matrix and get the reordering of the non-zero entries. More... | |
bool | is_transpose (const Sparsity &y) const |
Check if the sparsity is the transpose of another. More... | |
bool | is_reshape (const Sparsity &y) const |
Check if the sparsity is a reshape of another. More... | |
bool | is_subset (const Sparsity &rhs) const |
Is subset? More... | |
Sparsity | sparsity_cast_mod (const Sparsity &X, const Sparsity &Y) const |
Propagates subset according to sparsity cast. More... | |
Sparsity | pattern_inverse () const |
Take the inverse of a sparsity pattern; flip zeros and non-zeros. More... | |
void | enlarge (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &rr, const std::vector< casadi_int > &cc, bool ind1=false) |
Enlarge matrix. More... | |
void | enlargeRows (casadi_int nrow, const std::vector< casadi_int > &rr, bool ind1=false) |
Enlarge the matrix along the first dimension (i.e. insert rows) More... | |
void | enlargeColumns (casadi_int ncol, const std::vector< casadi_int > &cc, bool ind1=false) |
Enlarge the matrix along the second dimension (i.e. insert columns) More... | |
Sparsity | makeDense (std::vector< casadi_int > &mapping) const |
Make a patten dense. More... | |
std::vector< casadi_int > | erase (const std::vector< casadi_int > &rr, const std::vector< casadi_int > &cc, bool ind1=false) |
Erase rows and/or columns of a matrix. More... | |
std::vector< casadi_int > | erase (const std::vector< casadi_int > &rr, bool ind1=false) |
Erase elements of a matrix. More... | |
void | append (const Sparsity &sp) |
Append another sparsity patten vertically (NOTE: only efficient if vector) More... | |
void | appendColumns (const Sparsity &sp) |
Append another sparsity patten horizontally. More... | |
bool | is_scalar (bool scalar_and_dense=false) const |
Is scalar? More... | |
bool | is_dense () const |
Is dense? More... | |
bool | is_row () const |
Check if the pattern is a row vector (i.e. size1()==1) More... | |
bool | is_column () const |
Check if the pattern is a column vector (i.e. size2()==1) More... | |
bool | is_vector () const |
Check if the pattern is a row or column vector. More... | |
bool | is_diag () const |
Is diagonal? More... | |
bool | is_square () const |
Is square? More... | |
bool | is_symmetric () const |
Is symmetric? More... | |
bool | is_triu (bool strictly=false) const |
Is upper triangular? More... | |
bool | is_tril (bool strictly=false) const |
Is lower triangular? More... | |
bool | is_singular () const |
Check whether the sparsity-pattern indicates structural singularity. More... | |
bool | is_permutation () const |
Is this a permutation matrix? More... | |
bool | is_selection (bool allow_empty=false) const |
Is this a selection matrix? More... | |
bool | is_orthonormal (bool allow_empty=false) const |
Are both rows and columns orthonormal ? More... | |
bool | is_orthonormal_rows (bool allow_empty=false) const |
Are the rows of the pattern orthonormal ? More... | |
bool | is_orthonormal_columns (bool allow_empty=false) const |
Are the columns of the pattern orthonormal ? More... | |
bool | rowsSequential (bool strictly=true) const |
Do the rows appear sequentially on each column. More... | |
void | removeDuplicates (std::vector< casadi_int > &mapping) |
Remove duplicate entries. More... | |
std::vector< casadi_int > | etree (bool ata=false) const |
Calculate the elimination tree. More... | |
Sparsity | ldl (std::vector< casadi_int > &p, bool amd=true) const |
Symbolic LDL factorization. More... | |
void | qr_sparse (Sparsity &V, Sparsity &R, std::vector< casadi_int > &prinv, std::vector< casadi_int > &pc, bool amd=true) const |
Symbolic QR factorization. More... | |
casadi_int | dfs (casadi_int j, casadi_int top, std::vector< casadi_int > &xi, std::vector< casadi_int > &pstack, const std::vector< casadi_int > &pinv, std::vector< bool > &marked) const |
Depth-first search on the adjacency graph of the sparsity. More... | |
casadi_int | scc (std::vector< casadi_int > &index, std::vector< casadi_int > &offset) const |
Find the strongly connected components of the bigraph defined by the sparsity pattern. More... | |
casadi_int | btf (std::vector< casadi_int > &rowperm, std::vector< casadi_int > &colperm, std::vector< casadi_int > &rowblock, std::vector< casadi_int > &colblock, std::vector< casadi_int > &coarse_rowblock, std::vector< casadi_int > &coarse_colblock) const |
Calculate the block triangular form (BTF) More... | |
std::vector< casadi_int > | amd () const |
Approximate minimal degree preordering. More... | |
std::vector< casadi_int > | find (bool ind1=SWIG_IND1) const |
Get the location of all non-zero elements as they would appear in a Dense matrix. More... | |
Sparsity | uni_coloring (const Sparsity &AT=Sparsity(), casadi_int cutoff=std::numeric_limits< casadi_int >::max()) const |
Perform a unidirectional coloring: A greedy distance-2 coloring algorithm. More... | |
Sparsity | star_coloring (casadi_int ordering=1, casadi_int cutoff=std::numeric_limits< casadi_int >::max()) const |
Perform a star coloring of a symmetric matrix: More... | |
Sparsity | star_coloring2 (casadi_int ordering=1, casadi_int cutoff=std::numeric_limits< casadi_int >::max()) const |
Perform a star coloring of a symmetric matrix: More... | |
std::vector< casadi_int > | largest_first () const |
Order the columns by decreasing degree. More... | |
Sparsity | pmult (const std::vector< casadi_int > &p, bool permute_rows=true, bool permute_columns=true, bool invert_permutation=false) const |
Permute rows and/or columns. More... | |
std::string | dim (bool with_nz=false) const |
Get the dimension as a string. More... | |
std::string | postfix_dim () const |
Dimension string as a postfix to a name. More... | |
std::string | repr_el (casadi_int k) const |
Describe the nonzero location k as a string. More... | |
void | spy (std::ostream &stream=casadi::uout()) const |
Print a textual representation of sparsity. More... | |
void | spy_matlab (const std::string &mfile) const |
Generate a script for Matlab or Octave which visualizes. More... | |
void | export_code (const std::string &lang, std::ostream &stream=casadi::uout(), const Dict &options=Dict()) const |
Export matrix in specific language. More... | |
std::size_t | hash () const |
std::string | class_name () const |
Get class name. More... | |
void | disp (std::ostream &stream, bool more=false) const |
Print a description of the object. More... | |
std::string | get_str (bool more=false) const |
Get string representation. More... | |
bool | is_null () const |
Is a null pointer? More... | |
casadi_int | __hash__ () const |
Returns a number that is unique for a given Node. More... | |
Check if two sparsity patterns are identical | |
bool | is_equal (const Sparsity &y) const |
bool | is_equal (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &colind, const std::vector< casadi_int > &row) const |
bool | operator== (const Sparsity &y) const |
Size and element counting | |
casadi_int | size1 () const |
Get the number of rows. More... | |
casadi_int | rows () const |
Get the number of rows, Octave-style syntax. More... | |
casadi_int | size2 () const |
Get the number of columns. More... | |
casadi_int | columns () const |
Get the number of columns, Octave-style syntax. More... | |
casadi_int | numel () const |
The total number of elements, including structural zeros, i.e. size2()*size1() More... | |
double | density () const |
The percentage of nonzero. More... | |
bool | is_empty (bool both=false) const |
Check if the sparsity is empty. More... | |
casadi_int | nnz () const |
Get the number of (structural) non-zeros. More... | |
casadi_int | nnz_upper (bool strictly=false) const |
Number of non-zeros in the upper triangular half,. More... | |
casadi_int | nnz_lower (bool strictly=false) const |
Number of non-zeros in the lower triangular half,. More... | |
casadi_int | nnz_diag () const |
Number of non-zeros on the diagonal, i.e. the number of elements (i, j) with j==i. More... | |
casadi_int | bw_upper () const |
Upper half-bandwidth. More... | |
casadi_int | bw_lower () const |
Lower half-bandwidth. More... | |
std::pair< casadi_int, casadi_int > | size () const |
Get the shape. More... | |
casadi_int | size (casadi_int axis) const |
Get the size along a particular dimensions. More... | |
Sparsity | combine (const Sparsity &y, bool f0x_is_zero, bool function0_is_zero) const |
Combine two sparsity patterns. More... | |
Sparsity | unite (const Sparsity &y) const |
Union of two sparsity patterns. More... | |
Sparsity | operator+ (const Sparsity &b) const |
Union of two sparsity patterns. More... | |
Sparsity | intersect (const Sparsity &y) const |
Intersection of two sparsity patterns. More... | |
Sparsity | operator* (const Sparsity &b) const |
Intersection of two sparsity patterns. More... | |
Static Public Member Functions | |
static Sparsity | upper (casadi_int n) |
Create a upper triangular square sparsity pattern *. More... | |
static Sparsity | lower (casadi_int n) |
Create a lower triangular square sparsity pattern *. More... | |
static Sparsity | band (casadi_int n, casadi_int p) |
Create a single band in a square sparsity pattern. More... | |
static Sparsity | banded (casadi_int n, casadi_int p) |
Create banded square sparsity pattern. More... | |
static Sparsity | rowcol (const std::vector< casadi_int > &row, const std::vector< casadi_int > &col, casadi_int nrow, casadi_int ncol) |
Construct a block sparsity pattern from (row, col) vectors. More... | |
static Sparsity | nonzeros (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &nz, bool ind1=SWIG_IND1) |
Create a sparsity from nonzeros. More... | |
static Sparsity | permutation (const std::vector< casadi_int > &p, bool invert=false) |
Construct a permutation matrix P from a permutation vector p. More... | |
static Sparsity | from_file (const std::string &filename, const std::string &format_hint="") |
static Sparsity | deserialize (std::istream &stream) |
Build Sparsity from serialization. More... | |
static Sparsity | deserialize (const std::string &s) |
Build Sparsity from serialization. More... | |
static Sparsity | deserialize (DeserializingStream &s) |
Deserialize. More... | |
static std::string | type_name () |
Readable name of the public class. More... | |
static bool | test_cast (const SharedObjectInternal *ptr) |
Check if a particular cast is allowed. More... | |
static Sparsity | kkt (const Sparsity &H, const Sparsity &J, bool with_x_diag=true, bool with_lam_g_diag=true) |
Get KKT system sparsity. More... | |
static Sparsity | scalar (bool dense_scalar=true) |
Create a scalar sparsity pattern *. More... | |
static Sparsity | dense (casadi_int nrow, casadi_int ncol=1) |
Create a dense rectangular sparsity pattern *. More... | |
static Sparsity | dense (const std::pair< casadi_int, casadi_int > &rc) |
Create a dense rectangular sparsity pattern *. More... | |
static Sparsity | unit (casadi_int n, casadi_int el) |
Create the sparsity pattern for a unit vector of length n and a nonzero on. More... | |
static Sparsity | diag (casadi_int nrow) |
Create diagonal sparsity pattern *. More... | |
static Sparsity | diag (casadi_int nrow, casadi_int ncol) |
Create diagonal sparsity pattern *. More... | |
static Sparsity | diag (const std::pair< casadi_int, casadi_int > &rc) |
Create diagonal sparsity pattern *. More... | |
static Sparsity | triplet (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &row, const std::vector< casadi_int > &col, std::vector< casadi_int > &mapping, bool invert_mapping) |
Create a sparsity pattern given the nonzeros in sparse triplet form *. More... | |
static Sparsity | triplet (casadi_int nrow, casadi_int ncol, const std::vector< casadi_int > &row, const std::vector< casadi_int > &col) |
Create a sparsity pattern given the nonzeros in sparse triplet form *. More... | |
static Sparsity | compressed (const std::vector< casadi_int > &v, bool order_rows=false) |
|
explicit |
casadi::Sparsity::Sparsity | ( | casadi_int | nrow, |
casadi_int | ncol | ||
) |
Extra doc: https://github.com/casadi/casadi/wiki/L_ba
casadi::Sparsity::Sparsity | ( | casadi_int | nrow, |
casadi_int | ncol, | ||
const std::vector< casadi_int > & | colind, | ||
const std::vector< casadi_int > & | row, | ||
bool | order_rows = false |
||
) |
|
explicit |
Extra doc: https://github.com/casadi/casadi/wiki/L_bb
|
inherited |
If the Object does not point to any node, "0" is returned.
Extra doc: https://github.com/casadi/casadi/wiki/L_av
casadi_int casadi::Sparsity::add_nz | ( | casadi_int | rr, |
casadi_int | cc | ||
) |
Add the element if it does not exist and copy object if it's not unique
Extra doc: https://github.com/casadi/casadi/wiki/L_ce
std::vector<casadi_int> casadi::Sparsity::amd | ( | ) | const |
Fill-reducing ordering applied to the sparsity pattern of a linear system prior to factorization. The system must be symmetric, for an unsymmetric matrix A, first form the square of the pattern, A'*A.
The implementation is a modified version of cs_amd in CSparse Copyright(c) Timothy A. Davis, 2006-2009 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d8
void casadi::Sparsity::append | ( | const Sparsity & | sp | ) |
void casadi::Sparsity::appendColumns | ( | const Sparsity & | sp | ) |
|
static |
band(n, 0) is equivalent to diag(n)
band(n, -1) has a band below the diagonal
p | indicate |
Extra doc: https://github.com/casadi/casadi/wiki/L_bj
|
static |
banded(n, 0) is equivalent to diag(n)
banded(n, 1) is tri-diagonal matrix
Extra doc: https://github.com/casadi/casadi/wiki/L_bk
casadi_int casadi::Sparsity::btf | ( | std::vector< casadi_int > & | rowperm, |
std::vector< casadi_int > & | colperm, | ||
std::vector< casadi_int > & | rowblock, | ||
std::vector< casadi_int > & | colblock, | ||
std::vector< casadi_int > & | coarse_rowblock, | ||
std::vector< casadi_int > & | coarse_colblock | ||
) | const |
See Direct Methods for Sparse Linear Systems by Davis (2006).
The function computes the Dulmage-Mendelsohn decomposition, which allows you to reorder the rows and columns of a matrix to bring it into block triangular form (BTF).
It will not consider the distance of off-diagonal elements to the diagonal: there is no guarantee you will get a block-diagonal matrix if you supply a randomly permuted block-diagonal matrix.
If your matrix is symmetrical, this method is of limited use; permutation can make it non-symmetric.
The implementation is a modified version of cs_dmperm in CSparse Copyright(c) Timothy A. Davis, 2006-2009 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d7
casadi_int casadi::Sparsity::bw_lower | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_by
casadi_int casadi::Sparsity::bw_upper | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_bx
|
inherited |
Extra doc: https://github.com/casadi/casadi/wiki/L_au
casadi_int casadi::Sparsity::colind | ( | casadi_int | cc | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cb
|
inline |
Definition at line 333 of file sparsity.hpp.
Sparsity casadi::Sparsity::combine | ( | const Sparsity & | y, |
bool | f0x_is_zero, | ||
bool | function0_is_zero | ||
) | const |
Returns the new sparsity pattern as well as a mapping with the same length as the number of non-zero elements The mapping matrix contains the arguments for each nonzero, the first bit indicates if the first argument is nonzero, the second bit indicates if the second argument is nonzero (note that none of, one of or both of the arguments can be nonzero)
Extra doc: https://github.com/casadi/casadi/wiki/L_cl
std::vector<casadi_int> casadi::Sparsity::compress | ( | ) | const |
|
static |
Create from a single vector containing the pattern in compressed column storage format: The format: The first two entries are the number of rows (nrow) and columns (ncol) The next ncol+1 entries are the column offsets (colind). Note that the last element, colind[ncol], gives the number of nonzeros The last colind[ncol] entries are the row indices
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_be
|
inlinestatic |
Extra doc: https://github.com/casadi/casadi/wiki/L_be
Definition at line 155 of file sparsity.hpp.
double casadi::Sparsity::density | ( | ) | const |
Equivalent to (100.0 * nnz())/numel(), but avoids overflow
Extra doc: https://github.com/casadi/casadi/wiki/L_br
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_c4
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_c6
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_c3
casadi_int casadi::Sparsity::dfs | ( | casadi_int | j, |
casadi_int | top, | ||
std::vector< casadi_int > & | xi, | ||
std::vector< casadi_int > & | pstack, | ||
const std::vector< casadi_int > & | pinv, | ||
std::vector< bool > & | marked | ||
) | const |
See Direct Methods for Sparse Linear Systems by Davis (2006).
Extra doc: https://github.com/casadi/casadi/wiki/L_d5
|
inlinestatic |
Extra doc: https://github.com/casadi/casadi/wiki/L_bi
Definition at line 183 of file sparsity.hpp.
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bi
|
inlinestatic |
Extra doc: https://github.com/casadi/casadi/wiki/L_bi
Definition at line 185 of file sparsity.hpp.
std::string casadi::Sparsity::dim | ( | bool | with_nz = false | ) | const |
|
inherited |
void casadi::Sparsity::enlarge | ( | casadi_int | nrow, |
casadi_int | ncol, | ||
const std::vector< casadi_int > & | rr, | ||
const std::vector< casadi_int > & | cc, | ||
bool | ind1 = false |
||
) |
Make the matrix larger by inserting empty rows and columns, keeping the existing non-zeros
For the matrices A to B A(m, n) length(jj)=m , length(ii)=n B(nrow, ncol)
A=enlarge(m, n, ii, jj) makes sure that
B[jj, ii] == A
Extra doc: https://github.com/casadi/casadi/wiki/L_cr
void casadi::Sparsity::enlargeColumns | ( | casadi_int | ncol, |
const std::vector< casadi_int > & | cc, | ||
bool | ind1 = false |
||
) |
Extra doc: https://github.com/casadi/casadi/wiki/L_ct
void casadi::Sparsity::enlargeRows | ( | casadi_int | nrow, |
const std::vector< casadi_int > & | rr, | ||
bool | ind1 = false |
||
) |
Extra doc: https://github.com/casadi/casadi/wiki/L_cs
std::vector<casadi_int> casadi::Sparsity::erase | ( | const std::vector< casadi_int > & | rr, |
bool | ind1 = false |
||
) |
Extra doc: https://github.com/casadi/casadi/wiki/L_cw
std::vector<casadi_int> casadi::Sparsity::erase | ( | const std::vector< casadi_int > & | rr, |
const std::vector< casadi_int > & | cc, | ||
bool | ind1 = false |
||
) |
Extra doc: https://github.com/casadi/casadi/wiki/L_cv
std::vector<casadi_int> casadi::Sparsity::etree | ( | bool | ata = false | ) | const |
See Direct Methods for Sparse Linear Systems by Davis (2006). If the parameter ata is false, the algorithm is equivalent to MATLAB's etree(A), except that the indices are zero-based. If ata is true, the algorithm is equivalent to MATLAB's etree(A, 'col').
The implementation is a modified version of cs_etree in CSparse Copyright(c) Timothy A. Davis, 2006-2009 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d2
void casadi::Sparsity::export_code | ( | const std::string & | lang, |
std::ostream & | stream = casadi::uout() , |
||
const Dict & | options = Dict() |
||
) | const |
lang: only 'matlab' supported for now
* options: * inline: Indicates if you want everything on a single line (default: False) * name: Name of exported variable (default: 'sp') * as_matrix: Matlab does not have a sparsity object. (default: false) * With this option true, a numeric matrix will be constructed *
Extra doc: https://github.com/casadi/casadi/wiki/L_dj
std::vector<casadi_int> casadi::Sparsity::find | ( | bool | ind1 = SWIG_IND1 | ) | const |
A : DenseMatrix 4 x 3 B : SparseMatrix 4 x 3 , 5 structural non-zeros
k = A.find() A[k] will contain the elements of A that are non-zero in B
Inverse of nonzeros
.
Extra doc: https://github.com/casadi/casadi/wiki/L_da
|
static |
void casadi::Sparsity::get_ccs | ( | std::vector< casadi_int > & | colind, |
std::vector< casadi_int > & | row | ||
) | const |
std::vector<casadi_int> casadi::Sparsity::get_col | ( | ) | const |
Together with the row-vector, this vector gives the sparsity of the matrix in sparse triplet format, i.e. the column and row for each non-zero elements
Extra doc: https://github.com/casadi/casadi/wiki/L_cd
std::vector<casadi_int> casadi::Sparsity::get_colind | ( | ) | const |
Together with the row-vector, one obtains the sparsity pattern in the column compressed format.
Extra doc: https://github.com/casadi/casadi/wiki/L_ca
void casadi::Sparsity::get_crs | ( | std::vector< casadi_int > & | rowind, |
std::vector< casadi_int > & | col | ||
) | const |
Sparsity casadi::Sparsity::get_diag | ( | std::vector< casadi_int > & | mapping | ) | const |
Get the diagonal of the matrix/create a diagonal matrix (mapping will contain the nonzero mapping) When the input is square, the diagonal elements are returned. If the input is vector-like, a diagonal matrix is constructed with it.
std::vector<casadi_int> casadi::Sparsity::get_lower | ( | ) | const |
casadi_int casadi::Sparsity::get_nz | ( | casadi_int | rr, |
casadi_int | cc | ||
) | const |
return -1 if the element does not exist
Extra doc: https://github.com/casadi/casadi/wiki/L_cf
std::vector<casadi_int> casadi::Sparsity::get_nz | ( | const std::vector< casadi_int > & | rr, |
const std::vector< casadi_int > & | cc | ||
) | const |
return -1 if the element does not exist
Extra doc: https://github.com/casadi/casadi/wiki/L_cg
void casadi::Sparsity::get_nz | ( | std::vector< casadi_int > & | indices | ) | const |
The index vector is used both for input and outputs and must be sorted by increasing nonzero index, i.e. column-wise. Elements not found in the sparsity pattern are set to -1.
Extra doc: https://github.com/casadi/casadi/wiki/L_ch
std::vector<casadi_int> casadi::Sparsity::get_row | ( | ) | const |
Together with the column-vector, this vector gives the sparsity of the matrix in sparse triplet format, and together with the colind vector, one obtains the sparsity in column compressed format.
Extra doc: https://github.com/casadi/casadi/wiki/L_c9
|
inlineinherited |
Definition at line 138 of file shared_object.hpp.
void casadi::Sparsity::get_triplet | ( | std::vector< casadi_int > & | row, |
std::vector< casadi_int > & | col | ||
) | const |
std::vector<casadi_int> casadi::Sparsity::get_upper | ( | ) | const |
bool casadi::Sparsity::has_nz | ( | casadi_int | rr, |
casadi_int | cc | ||
) | const |
std::size_t casadi::Sparsity::hash | ( | ) | const |
Dict casadi::Sparsity::info | ( | ) | const |
Obtain information about sparsity
Returns the new sparsity pattern as well as a mapping with the same length as the number of non-zero elements The value is 1 if the non-zero comes from the first (i.e. this) object, 2 if it is from the second and 3 (i.e. 1 | 2) if from both
Extra doc: https://github.com/casadi/casadi/wiki/L_cn
bool casadi::Sparsity::is_column | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cy
bool casadi::Sparsity::is_dense | ( | ) | const |
bool casadi::Sparsity::is_diag | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_247
bool casadi::Sparsity::is_empty | ( | bool | both = false | ) | const |
A sparsity is considered empty if one of the dimensions is zero (or optionally both dimensions)
Extra doc: https://github.com/casadi/casadi/wiki/L_bs
bool casadi::Sparsity::is_equal | ( | casadi_int | nrow, |
casadi_int | ncol, | ||
const std::vector< casadi_int > & | colind, | ||
const std::vector< casadi_int > & | row | ||
) | const |
bool casadi::Sparsity::is_equal | ( | const Sparsity & | y | ) | const |
|
inherited |
bool casadi::Sparsity::is_orthonormal | ( | bool | allow_empty = false | ) | const |
[in] | allow_empty | Disregard empty rows and columns in the analysis |
Extra doc: https://github.com/casadi/casadi/wiki/L_24f
bool casadi::Sparsity::is_orthonormal_columns | ( | bool | allow_empty = false | ) | const |
[in] | allow_empty | Disregard empty columns in the analysis |
Extra doc: https://github.com/casadi/casadi/wiki/L_24h
bool casadi::Sparsity::is_orthonormal_rows | ( | bool | allow_empty = false | ) | const |
[in] | allow_empty | Disregard empty rows in the analysis |
Extra doc: https://github.com/casadi/casadi/wiki/L_24g
bool casadi::Sparsity::is_permutation | ( | ) | const |
A Matrix P is permutation matrix if right multiplication with a dense vector v leads to a vector with the same elements, but permuted.
Implies square
Equivalent to is_orthonormal(false)
Extra doc: https://github.com/casadi/casadi/wiki/L_24d
bool casadi::Sparsity::is_reshape | ( | const Sparsity & | y | ) | const |
bool casadi::Sparsity::is_row | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cx
bool casadi::Sparsity::is_scalar | ( | bool | scalar_and_dense = false | ) | const |
bool casadi::Sparsity::is_selection | ( | bool | allow_empty = false | ) | const |
A Matrix S is selection matrix if right multiplication with a dense vector leads to a vector with a subset of the elements of the original vector
[in] | allow_empty | Allow the resultant vector to have structural zeros |
Equivalent to is_orthonormal_rows(allow_empty)
Extra doc: https://github.com/casadi/casadi/wiki/L_24e
bool casadi::Sparsity::is_singular | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_24c
bool casadi::Sparsity::is_square | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_248
bool casadi::Sparsity::is_stacked | ( | const Sparsity & | y, |
casadi_int | n | ||
) | const |
bool casadi::Sparsity::is_subset | ( | const Sparsity & | rhs | ) | const |
bool casadi::Sparsity::is_symmetric | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_249
bool casadi::Sparsity::is_transpose | ( | const Sparsity & | y | ) | const |
bool casadi::Sparsity::is_tril | ( | bool | strictly = false | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_24b
bool casadi::Sparsity::is_triu | ( | bool | strictly = false | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_24a
bool casadi::Sparsity::is_vector | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cz
|
static |
[H + I1, J'; J, I2] where I1 and I2 are optional
Extra doc: https://github.com/casadi/casadi/wiki/L_dk
std::vector<casadi_int> casadi::Sparsity::largest_first | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_de
Sparsity casadi::Sparsity::ldl | ( | std::vector< casadi_int > & | p, |
bool | amd = true |
||
) | const |
Returns the sparsity pattern of L^T
The implementation is a modified version of LDL Copyright(c) Timothy A. Davis, 2005-2013 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d3
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bh
Sparsity casadi::Sparsity::makeDense | ( | std::vector< casadi_int > & | mapping | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cu
casadi_int casadi::Sparsity::nnz | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_bt
casadi_int casadi::Sparsity::nnz_diag | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_bw
casadi_int casadi::Sparsity::nnz_lower | ( | bool | strictly = false | ) | const |
i.e. the number of elements (i, j) with j<=i
Extra doc: https://github.com/casadi/casadi/wiki/L_bv
casadi_int casadi::Sparsity::nnz_upper | ( | bool | strictly = false | ) | const |
i.e. the number of elements (i, j) with j>=i
Extra doc: https://github.com/casadi/casadi/wiki/L_bu
|
static |
Inverse of find()
Extra doc: https://github.com/casadi/casadi/wiki/L_23g
casadi_int casadi::Sparsity::numel | ( | ) | const |
|
inline |
Definition at line 295 of file sparsity.hpp.
Returns the new sparsity pattern as well as a mapping with the same length as the number of non-zero elements The value is 1 if the non-zero comes from the first (i.e. this) object, 2 if it is from the second and 3 (i.e. 1 | 2) if from both
Extra doc: https://github.com/casadi/casadi/wiki/L_cn
Extra doc: https://github.com/casadi/casadi/wiki/L_cm
|
inline |
Definition at line 291 of file sparsity.hpp.
Sparsity casadi::Sparsity::pattern_inverse | ( | ) | const |
|
static |
Right multiplication of P with a vector leads to the same results as indexing that vector with p
P @ v = v[p]
The inverse of a permutation matrix is equal to its transpose (property of orthonormality)
Extra doc: https://github.com/casadi/casadi/wiki/L_244
const std::vector<casadi_int> casadi::Sparsity::permutation_vector | ( | bool | invert = false | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_245
Sparsity casadi::Sparsity::pmult | ( | const std::vector< casadi_int > & | p, |
bool | permute_rows = true , |
||
bool | permute_columns = true , |
||
bool | invert_permutation = false |
||
) | const |
Multiply the sparsity with a permutation matrix from the left and/or from the right P * A * trans(P), A * trans(P) or A * trans(P) with P defined by an index vector containing the row for each col. As an alternative, P can be transposed (inverted).
Extra doc: https://github.com/casadi/casadi/wiki/L_df
std::string casadi::Sparsity::postfix_dim | ( | ) | const |
Rules:
Otherwise: "[5x10,3nz]"
Extra doc: https://github.com/casadi/casadi/wiki/L_dg
void casadi::Sparsity::qr_sparse | ( | Sparsity & | V, |
Sparsity & | R, | ||
std::vector< casadi_int > & | prinv, | ||
std::vector< casadi_int > & | pc, | ||
bool | amd = true |
||
) | const |
Returns the sparsity pattern of V (compact representation of Q) and R as well as vectors needed for the numerical factorization and solution. The implementation is a modified version of CSparse Copyright(c) Timothy A. Davis, 2006-2009 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d4
void casadi::Sparsity::removeDuplicates | ( | std::vector< casadi_int > & | mapping | ) |
The same indices will be removed from the mapping vector, which must have the same length as the number of nonzeros
Extra doc: https://github.com/casadi/casadi/wiki/L_d1
std::string casadi::Sparsity::repr_el | ( | casadi_int | k | ) | const |
void casadi::Sparsity::resize | ( | casadi_int | nrow, |
casadi_int | ncol | ||
) |
casadi_int casadi::Sparsity::row | ( | casadi_int | el | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_cc
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bl
|
inline |
Definition at line 327 of file sparsity.hpp.
bool casadi::Sparsity::rowsSequential | ( | bool | strictly = true | ) | const |
[in] | strictly | if true, then do not allow multiple entries |
Extra doc: https://github.com/casadi/casadi/wiki/L_d0
|
inlinestatic |
Extra doc: https://github.com/casadi/casadi/wiki/L_bd
Definition at line 146 of file sparsity.hpp.
casadi_int casadi::Sparsity::scc | ( | std::vector< casadi_int > & | index, |
std::vector< casadi_int > & | offset | ||
) | const |
of a square matrix
See Direct Methods for Sparse Linear Systems by Davis (2006). Returns:
In the case that the matrix is symmetric, the result has a particular interpretation: Given a symmetric matrix A and n = A.scc(p, r)
=> A[p, p] will appear block-diagonal with n blocks and with the indices of the block boundaries to be found in r.
The implementation is a modified version of cs_scc in CSparse Copyright(c) Timothy A. Davis, 2006-2009 Licensed as a derivative work under the GNU LGPL
Extra doc: https://github.com/casadi/casadi/wiki/L_d6
std::string casadi::Sparsity::serialize | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_c2
void casadi::Sparsity::serialize | ( | SerializingStream & | s | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_c5
std::pair<casadi_int, casadi_int> casadi::Sparsity::size | ( | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_bz
casadi_int casadi::Sparsity::size | ( | casadi_int | axis | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_c0
casadi_int casadi::Sparsity::size1 | ( | ) | const |
casadi_int casadi::Sparsity::size2 | ( | ) | const |
Assumption: 'this' is a subset of X.
Example:
X = [ * . . .; . * . .; . . * .; . . . *] Y = [* *; . . ; * *] x = [ * . . .; . . . .; . . * .; . . . *] returns [* *; . . ; . *]
Extra doc: https://github.com/casadi/casadi/wiki/L_246
void casadi::Sparsity::spy | ( | std::ostream & | stream = casadi::uout() | ) | const |
Extra doc: https://github.com/casadi/casadi/wiki/L_dh
void casadi::Sparsity::spy_matlab | ( | const std::string & | mfile | ) | const |
the sparsity using the spy command
Extra doc: https://github.com/casadi/casadi/wiki/L_di
Sparsity casadi::Sparsity::star_coloring | ( | casadi_int | ordering = 1 , |
casadi_int | cutoff = std::numeric_limits< casadi_int >::max() |
||
) | const |
A greedy distance-2 coloring algorithm Algorithm 4.1 in What Color Is Your Jacobian? Graph Coloring for Computing Derivatives A. H. GEBREMEDHIN, F. MANNE, A. POTHEN SIAM Rev., 47(4), 629–705 (2006)
Ordering options: None (0), largest first (1)
Extra doc: https://github.com/casadi/casadi/wiki/L_dc
Sparsity casadi::Sparsity::star_coloring2 | ( | casadi_int | ordering = 1 , |
casadi_int | cutoff = std::numeric_limits< casadi_int >::max() |
||
) | const |
A new greedy distance-2 coloring algorithm Algorithm 4.1 in NEW ACYCLIC AND STAR COLORING ALGORITHMS WITH APPLICATION TO COMPUTING HESSIANS A. H. GEBREMEDHIN, A. TARAFDAR, F. MANNE, A. POTHEN SIAM J. SCI. COMPUT. Vol. 29, No. 3, pp. 1042–1072 (2007)
Ordering options: None (0), largest first (1)
Extra doc: https://github.com/casadi/casadi/wiki/L_dd
Sparsity casadi::Sparsity::sub | ( | const std::vector< casadi_int > & | rr, |
const Sparsity & | sp, | ||
std::vector< casadi_int > & | mapping, | ||
bool | ind1 = false |
||
) | const |
Returns the sparsity of the corresponding elements, with a mapping such that submatrix[k] = originalmatrix[mapping[k]]
Extra doc: https://github.com/casadi/casadi/wiki/L_cj
Sparsity casadi::Sparsity::sub | ( | const std::vector< casadi_int > & | rr, |
const std::vector< casadi_int > & | cc, | ||
std::vector< casadi_int > & | mapping, | ||
bool | ind1 = false |
||
) | const |
Returns the sparsity of the submatrix, with a mapping such that submatrix[k] = originalmatrix[mapping[k]]
Extra doc: https://github.com/casadi/casadi/wiki/L_ci
Sparsity casadi::Sparsity::T | ( | ) | const |
|
static |
void casadi::Sparsity::to_file | ( | const std::string & | filename, |
const std::string & | format_hint = "" |
||
) | const |
Sparsity casadi::Sparsity::transpose | ( | std::vector< casadi_int > & | mapping, |
bool | invert_mapping = false |
||
) | const |
[out] | mapping | the non-zeros of the original matrix for each non-zero of the new matrix |
Extra doc: https://github.com/casadi/casadi/wiki/L_ck
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bm
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bm
|
inlinestatic |
Definition at line 1128 of file sparsity.hpp.
Sparsity casadi::Sparsity::uni_coloring | ( | const Sparsity & | AT = Sparsity() , |
casadi_int | cutoff = std::numeric_limits< casadi_int >::max() |
||
) | const |
(Algorithm 3.1 in A. H. GEBREMEDHIN, F. MANNE, A. POTHEN)
Extra doc: https://github.com/casadi/casadi/wiki/L_db
|
static |
position el *
Extra doc: https://github.com/casadi/casadi/wiki/L_bf
Extra doc: https://github.com/casadi/casadi/wiki/L_cm
|
static |
Extra doc: https://github.com/casadi/casadi/wiki/L_bg