List of all members | Public Member Functions | Static Public Member Functions
casadi::Sparsity Class Reference

General sparsity class. More...

#include <sparsity.hpp>

Detailed Description

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:

  1. "colind" [length size2()+1], which contains the index to the first non-zero element on or after the corresponding column. All the non-zero elements of a particular column i are thus the elements with index el that fulfills: colind[i] <= el < colind[i+1].
  2. "row" [same length as the number of non-zero elements, nnz()] The rows for each of the structural non-zeros.

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

See also
Matrix
Author
Joel Andersson
Date
2010-2015

Extra doc: https://github.com/casadi/casadi/wiki/L_b9

Definition at line 96 of file sparsity.hpp.

Inheritance diagram for casadi::Sparsity:
Inheritance graph
[legend]

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)
 

Constructor & Destructor Documentation

◆ Sparsity() [1/4]

casadi::Sparsity::Sparsity ( casadi_int  dummy = 0)
explicit

◆ Sparsity() [2/4]

casadi::Sparsity::Sparsity ( casadi_int  nrow,
casadi_int  ncol 
)

◆ Sparsity() [3/4]

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 
)

◆ Sparsity() [4/4]

casadi::Sparsity::Sparsity ( const std::pair< casadi_int, casadi_int > &  rc)
explicit

Member Function Documentation

◆ __hash__()

casadi_int casadi::SharedObject::__hash__ ( ) const
inherited

If the Object does not point to any node, "0" is returned.

Extra doc: https://github.com/casadi/casadi/wiki/L_av

◆ add_nz()

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

◆ amd()

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

◆ append()

void casadi::Sparsity::append ( const Sparsity sp)

◆ appendColumns()

void casadi::Sparsity::appendColumns ( const Sparsity sp)

◆ band()

static Sparsity casadi::Sparsity::band ( casadi_int  n,
casadi_int  p 
)
static

band(n, 0) is equivalent to diag(n)
band(n, -1) has a band below the diagonal

Parameters
pindicate

Extra doc: https://github.com/casadi/casadi/wiki/L_bj

◆ banded()

static Sparsity casadi::Sparsity::banded ( casadi_int  n,
casadi_int  p 
)
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

◆ btf()

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.

See also
scc

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

Examples
matrix/btf.py.

◆ bw_lower()

casadi_int casadi::Sparsity::bw_lower ( ) const

◆ bw_upper()

casadi_int casadi::Sparsity::bw_upper ( ) const

◆ class_name()

std::string casadi::SharedObject::class_name ( ) const
inherited

◆ colind()

casadi_int casadi::Sparsity::colind ( casadi_int  cc) const

◆ columns()

casadi_int casadi::Sparsity::columns ( ) const
inline

Definition at line 333 of file sparsity.hpp.

◆ combine()

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

◆ compress()

std::vector<casadi_int> casadi::Sparsity::compress ( ) const

◆ compressed()

static Sparsity casadi::Sparsity::compressed ( const std::vector< casadi_int > &  v,
bool  order_rows = false 
)
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

◆ dense() [1/2]

static Sparsity casadi::Sparsity::dense ( casadi_int  nrow,
casadi_int  ncol = 1 
)
static

◆ dense() [2/2]

static Sparsity casadi::Sparsity::dense ( const std::pair< casadi_int, casadi_int > &  rc)
inlinestatic

Extra doc: https://github.com/casadi/casadi/wiki/L_be

Definition at line 155 of file sparsity.hpp.

◆ density()

double casadi::Sparsity::density ( ) const

Equivalent to (100.0 * nnz())/numel(), but avoids overflow

Extra doc: https://github.com/casadi/casadi/wiki/L_br

◆ deserialize() [1/3]

static Sparsity casadi::Sparsity::deserialize ( const std::string &  s)
static

◆ deserialize() [2/3]

static Sparsity casadi::Sparsity::deserialize ( DeserializingStream s)
static

◆ deserialize() [3/3]

static Sparsity casadi::Sparsity::deserialize ( std::istream &  stream)
static

◆ dfs()

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

◆ diag() [1/3]

static Sparsity casadi::Sparsity::diag ( casadi_int  nrow)
inlinestatic

Extra doc: https://github.com/casadi/casadi/wiki/L_bi

Definition at line 183 of file sparsity.hpp.

◆ diag() [2/3]

static Sparsity casadi::Sparsity::diag ( casadi_int  nrow,
casadi_int  ncol 
)
static

◆ diag() [3/3]

static Sparsity casadi::Sparsity::diag ( const std::pair< casadi_int, casadi_int > &  rc)
inlinestatic

Extra doc: https://github.com/casadi/casadi/wiki/L_bi

Definition at line 185 of file sparsity.hpp.

◆ dim()

std::string casadi::Sparsity::dim ( bool  with_nz = false) const

◆ disp()

void casadi::SharedObject::disp ( std::ostream &  stream,
bool  more = false 
) const
inherited

◆ enlarge()

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

◆ enlargeColumns()

void casadi::Sparsity::enlargeColumns ( casadi_int  ncol,
const std::vector< casadi_int > &  cc,
bool  ind1 = false 
)

◆ enlargeRows()

void casadi::Sparsity::enlargeRows ( casadi_int  nrow,
const std::vector< casadi_int > &  rr,
bool  ind1 = false 
)

◆ erase() [1/2]

std::vector<casadi_int> casadi::Sparsity::erase ( const std::vector< casadi_int > &  rr,
bool  ind1 = false 
)

◆ erase() [2/2]

std::vector<casadi_int> casadi::Sparsity::erase ( const std::vector< casadi_int > &  rr,
const std::vector< casadi_int > &  cc,
bool  ind1 = false 
)

◆ etree()

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

◆ export_code()

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

◆ find()

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

◆ from_file()

static Sparsity casadi::Sparsity::from_file ( const std::string &  filename,
const std::string &  format_hint = "" 
)
static

◆ get_ccs()

void casadi::Sparsity::get_ccs ( std::vector< casadi_int > &  colind,
std::vector< casadi_int > &  row 
) const

◆ get_col()

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

◆ get_colind()

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

◆ get_crs()

void casadi::Sparsity::get_crs ( std::vector< casadi_int > &  rowind,
std::vector< casadi_int > &  col 
) const

◆ get_diag()

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.

◆ get_lower()

std::vector<casadi_int> casadi::Sparsity::get_lower ( ) const

◆ get_nz() [1/3]

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

◆ get_nz() [2/3]

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

◆ get_nz() [3/3]

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

◆ get_row()

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

◆ get_str()

std::string casadi::SharedObject::get_str ( bool  more = false) const
inlineinherited

Definition at line 138 of file shared_object.hpp.

◆ get_triplet()

void casadi::Sparsity::get_triplet ( std::vector< casadi_int > &  row,
std::vector< casadi_int > &  col 
) const

◆ get_upper()

std::vector<casadi_int> casadi::Sparsity::get_upper ( ) const

◆ has_nz()

bool casadi::Sparsity::has_nz ( casadi_int  rr,
casadi_int  cc 
) const

◆ hash()

std::size_t casadi::Sparsity::hash ( ) const

◆ info()

Dict casadi::Sparsity::info ( ) const

Obtain information about sparsity

◆ intersect()

Sparsity casadi::Sparsity::intersect ( const Sparsity y) const

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

◆ is_column()

bool casadi::Sparsity::is_column ( ) const

◆ is_dense()

bool casadi::Sparsity::is_dense ( ) const

◆ is_diag()

bool casadi::Sparsity::is_diag ( ) const

◆ is_empty()

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

◆ is_equal() [1/2]

bool casadi::Sparsity::is_equal ( casadi_int  nrow,
casadi_int  ncol,
const std::vector< casadi_int > &  colind,
const std::vector< casadi_int > &  row 
) const

◆ is_equal() [2/2]

bool casadi::Sparsity::is_equal ( const Sparsity y) const

◆ is_null()

bool casadi::SharedObject::is_null ( ) const
inherited

◆ is_orthonormal()

bool casadi::Sparsity::is_orthonormal ( bool  allow_empty = false) const
Parameters
[in]allow_emptyDisregard empty rows and columns in the analysis

Extra doc: https://github.com/casadi/casadi/wiki/L_24f

◆ is_orthonormal_columns()

bool casadi::Sparsity::is_orthonormal_columns ( bool  allow_empty = false) const
Parameters
[in]allow_emptyDisregard empty columns in the analysis

Extra doc: https://github.com/casadi/casadi/wiki/L_24h

◆ is_orthonormal_rows()

bool casadi::Sparsity::is_orthonormal_rows ( bool  allow_empty = false) const
Parameters
[in]allow_emptyDisregard empty rows in the analysis

Extra doc: https://github.com/casadi/casadi/wiki/L_24g

◆ is_permutation()

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

◆ is_reshape()

bool casadi::Sparsity::is_reshape ( const Sparsity y) const

◆ is_row()

bool casadi::Sparsity::is_row ( ) const

◆ is_scalar()

bool casadi::Sparsity::is_scalar ( bool  scalar_and_dense = false) const

◆ is_selection()

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

Parameters
[in]allow_emptyAllow the resultant vector to have structural zeros

Equivalent to is_orthonormal_rows(allow_empty)

Extra doc: https://github.com/casadi/casadi/wiki/L_24e

◆ is_singular()

bool casadi::Sparsity::is_singular ( ) const

◆ is_square()

bool casadi::Sparsity::is_square ( ) const

◆ is_stacked()

bool casadi::Sparsity::is_stacked ( const Sparsity y,
casadi_int  n 
) const

◆ is_subset()

bool casadi::Sparsity::is_subset ( const Sparsity rhs) const

◆ is_symmetric()

bool casadi::Sparsity::is_symmetric ( ) const

◆ is_transpose()

bool casadi::Sparsity::is_transpose ( const Sparsity y) const

◆ is_tril()

bool casadi::Sparsity::is_tril ( bool  strictly = false) const

◆ is_triu()

bool casadi::Sparsity::is_triu ( bool  strictly = false) const

◆ is_vector()

bool casadi::Sparsity::is_vector ( ) const

◆ kkt()

static Sparsity casadi::Sparsity::kkt ( const Sparsity H,
const Sparsity J,
bool  with_x_diag = true,
bool  with_lam_g_diag = true 
)
static

[H + I1, J'; J, I2] where I1 and I2 are optional

Extra doc: https://github.com/casadi/casadi/wiki/L_dk

◆ largest_first()

std::vector<casadi_int> casadi::Sparsity::largest_first ( ) const

◆ ldl()

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

◆ lower()

static Sparsity casadi::Sparsity::lower ( casadi_int  n)
static

◆ makeDense()

Sparsity casadi::Sparsity::makeDense ( std::vector< casadi_int > &  mapping) const

◆ nnz()

casadi_int casadi::Sparsity::nnz ( ) const

◆ nnz_diag()

casadi_int casadi::Sparsity::nnz_diag ( ) const

◆ nnz_lower()

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

◆ nnz_upper()

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

◆ nonzeros()

static Sparsity casadi::Sparsity::nonzeros ( casadi_int  nrow,
casadi_int  ncol,
const std::vector< casadi_int > &  nz,
bool  ind1 = SWIG_IND1 
)
static

◆ numel()

casadi_int casadi::Sparsity::numel ( ) const

Beware of overflow

See also
nnz()

Extra doc: https://github.com/casadi/casadi/wiki/L_bq

◆ operator!=()

bool casadi::Sparsity::operator!= ( const Sparsity y) const
inline

Definition at line 295 of file sparsity.hpp.

◆ operator*()

Sparsity casadi::Sparsity::operator* ( const Sparsity b) const

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

◆ operator+()

Sparsity casadi::Sparsity::operator+ ( const Sparsity b) const

◆ operator==()

bool casadi::Sparsity::operator== ( const Sparsity y) const
inline

Definition at line 291 of file sparsity.hpp.

◆ pattern_inverse()

Sparsity casadi::Sparsity::pattern_inverse ( ) const

◆ permutation()

static Sparsity casadi::Sparsity::permutation ( const std::vector< casadi_int > &  p,
bool  invert = false 
)
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

◆ permutation_vector()

const std::vector<casadi_int> casadi::Sparsity::permutation_vector ( bool  invert = false) const

◆ pmult()

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

◆ postfix_dim()

std::string casadi::Sparsity::postfix_dim ( ) const

Rules:

  1. Dense and scalar: ""
  2. 0-by-0: "[]"
  3. Dense column vector: "[5]"
  4. Dense matrix: "[5x10]"
  5. Otherwise: "[5x10,3nz]"

    Extra doc: https://github.com/casadi/casadi/wiki/L_dg

◆ qr_sparse()

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

◆ removeDuplicates()

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

◆ repr_el()

std::string casadi::Sparsity::repr_el ( casadi_int  k) const

◆ resize()

void casadi::Sparsity::resize ( casadi_int  nrow,
casadi_int  ncol 
)

◆ row()

casadi_int casadi::Sparsity::row ( casadi_int  el) const

◆ rowcol()

static Sparsity casadi::Sparsity::rowcol ( const std::vector< casadi_int > &  row,
const std::vector< casadi_int > &  col,
casadi_int  nrow,
casadi_int  ncol 
)
static

◆ rows()

casadi_int casadi::Sparsity::rows ( ) const
inline

Definition at line 327 of file sparsity.hpp.

◆ rowsSequential()

bool casadi::Sparsity::rowsSequential ( bool  strictly = true) const
Parameters
[in]strictlyif true, then do not allow multiple entries

Extra doc: https://github.com/casadi/casadi/wiki/L_d0

◆ scalar()

static Sparsity casadi::Sparsity::scalar ( bool  dense_scalar = true)
inlinestatic

Extra doc: https://github.com/casadi/casadi/wiki/L_bd

Definition at line 146 of file sparsity.hpp.

◆ scc()

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:

  • Number of components
  • Offset for each components (length: 1 + number of components)
  • Indices for each components, component i has indices index[offset[i]], ..., index[offset[i+1]]

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

Examples
matrix/scc.py.

◆ serialize() [1/2]

std::string casadi::Sparsity::serialize ( ) const

◆ serialize() [2/2]

void casadi::Sparsity::serialize ( SerializingStream s) const

◆ size() [1/2]

std::pair<casadi_int, casadi_int> casadi::Sparsity::size ( ) const

◆ size() [2/2]

casadi_int casadi::Sparsity::size ( casadi_int  axis) const

◆ size1()

casadi_int casadi::Sparsity::size1 ( ) const

◆ size2()

casadi_int casadi::Sparsity::size2 ( ) const

◆ sparsity_cast_mod()

Sparsity casadi::Sparsity::sparsity_cast_mod ( const Sparsity X,
const Sparsity Y 
) const

Assumption: 'this' is a subset of X.

Example:

X = [ * . . .; . * . .; . . * .; . . . *] Y = [* *; . . ; * *] x = [ * . . .; . . . .; . . * .; . . . *] returns [* *; . . ; . *]

Extra doc: https://github.com/casadi/casadi/wiki/L_246

◆ spy()

void casadi::Sparsity::spy ( std::ostream &  stream = casadi::uout()) const

◆ spy_matlab()

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

◆ star_coloring()

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

◆ star_coloring2()

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

◆ sub() [1/2]

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

◆ sub() [2/2]

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

◆ T()

Sparsity casadi::Sparsity::T ( ) const

◆ test_cast()

static bool casadi::Sparsity::test_cast ( const SharedObjectInternal *  ptr)
static

◆ to_file()

void casadi::Sparsity::to_file ( const std::string &  filename,
const std::string &  format_hint = "" 
) const

Export sparsity pattern to file

Supported formats:

◆ transpose()

Sparsity casadi::Sparsity::transpose ( std::vector< casadi_int > &  mapping,
bool  invert_mapping = false 
) const
Parameters
[out]mappingthe non-zeros of the original matrix for each non-zero of the new matrix

Extra doc: https://github.com/casadi/casadi/wiki/L_ck

◆ triplet() [1/2]

static Sparsity casadi::Sparsity::triplet ( casadi_int  nrow,
casadi_int  ncol,
const std::vector< casadi_int > &  row,
const std::vector< casadi_int > &  col 
)
static

◆ triplet() [2/2]

static Sparsity casadi::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 
)
static

◆ type_name()

static std::string casadi::Sparsity::type_name ( )
inlinestatic

Definition at line 1128 of file sparsity.hpp.

◆ uni_coloring()

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

Examples
matrix/coloring.py.

◆ unit()

static Sparsity casadi::Sparsity::unit ( casadi_int  n,
casadi_int  el 
)
static

◆ unite()

Sparsity casadi::Sparsity::unite ( const Sparsity y) const

◆ upper()

static Sparsity casadi::Sparsity::upper ( casadi_int  n)
static

The documentation for this class was generated from the following file: