List of all members | Public Member Functions | Static Public Member Functions | Friends
casadi::GenericMatrix< MatType > Class Template Reference

Matrix base class. More...

#include <generic_matrix.hpp>

Detailed Description

template<typename MatType>
class casadi::GenericMatrix< MatType >

This is a common base class for MX and Matrix<>, introducing a uniform syntax and implementing common functionality using the curiously recurring template pattern (CRTP) idiom.
The class is designed with the idea that "everything is a matrix", that is, also scalars and vectors.
This philosophy makes it easy to use and to interface in particularly with Python and Matlab/Octave.
The syntax tries to stay as close as possible to the ublas syntax when it comes to vector/matrix operations.
Index starts with 0.
Index vec happens as follows: (rr, cc) -> k = rr+cc*size1()
Vectors are column vectors.
The storage format is Compressed Column Storage (CCS), similar to that used for sparse matrices in Matlab,
but unlike this format, we do allow for elements to be structurally non-zero but numerically zero.
The sparsity pattern, which is reference counted and cached, can be accessed with Sparsity& sparsity()

Author
Joel Andersson
Date
2012

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

Definition at line 75 of file generic_matrix.hpp.

Inheritance diagram for casadi::GenericMatrix< MatType >:
Inheritance graph
[legend]
Collaboration diagram for casadi::GenericMatrix< MatType >:
Collaboration graph
[legend]

Public Member Functions

casadi_int nnz () const
 Get the number of (structural) non-zero elements. More...
 
casadi_int nnz_lower () const
 Get the number of non-zeros in the lower triangular half. More...
 
casadi_int nnz_upper () const
 Get the number of non-zeros in the upper triangular half. More...
 
casadi_int nnz_diag () const
 Get get the number of non-zeros on the diagonal. More...
 
casadi_int numel () const
 Get the number of elements. More...
 
casadi_int size1 () const
 Get the first dimension (i.e. number of rows) More...
 
casadi_int rows () const
 Get the number of rows, Octave-style syntax. More...
 
casadi_int size2 () const
 Get the second dimension (i.e. number of columns) More...
 
casadi_int columns () const
 Get the number of columns, Octave-style syntax. More...
 
std::string dim (bool with_nz=false) const
 Get string representation of dimensions. 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...
 
bool is_empty (bool both=false) const
 Check if the sparsity is empty, i.e. if one of the dimensions is zero. More...
 
bool is_dense () const
 Check if the matrix expression is dense. More...
 
bool is_scalar (bool scalar_and_dense=false) const
 Check if the matrix expression is scalar. More...
 
bool is_square () const
 Check if the matrix expression is square. More...
 
bool is_vector () const
 Check if the matrix is a row or column vector. More...
 
bool is_row () const
 Check if the matrix is a row vector (i.e. size1()==1) More...
 
bool is_column () const
 Check if the matrix is a column vector (i.e. size2()==1) More...
 
bool is_triu () const
 Check if the matrix is upper triangular. More...
 
bool is_tril () const
 Check if the matrix is lower triangular. More...
 
Sparsity sparsity () const
 Get the sparsity pattern. More...
 
template<typename K >
const MatType nz (const K &k) const
 Get vector nonzero or slice of nonzeros. More...
 
template<typename K >
NonZeros< MatType, K > nz (const K &k)
 Access vector nonzero or slice of nonzeros. More...
 
template<typename RR >
const MatType operator() (const RR &rr) const
 Get vector element or slice. More...
 
template<typename RR , typename CC >
const MatType operator() (const RR &rr, const CC &cc) const
 Get Matrix element or slice. More...
 
template<typename RR >
SubIndex< MatType, RR > operator() (const RR &rr)
 Access Matrix elements (one argument) More...
 
template<typename RR , typename CC >
SubMatrix< MatType, RR, CC > operator() (const RR &rr, const CC &cc)
 Access Matrix elements (two arguments) More...
 
std::vector< casadi_int > get_row () const
 Get the sparsity pattern. See the Sparsity class for details. More...
 
std::vector< casadi_int > get_colind () const
 Get the sparsity pattern. See the Sparsity class for details. More...
 
const casadi_int * row () const
 Get the sparsity pattern. See the Sparsity class for details. More...
 
const casadi_int * colind () const
 Get the sparsity pattern. See the Sparsity class for details. More...
 
casadi_int row (casadi_int el) const
 Get the sparsity pattern. See the Sparsity class for details. More...
 
casadi_int colind (casadi_int col) const
 Get the sparsity pattern. See the Sparsity class for details. More...
 

Static Public Member Functions

static MatType logsumexp (const MatType &x)
 
static MatType interp1d (const std::vector< double > &x, const MatType &v, const std::vector< double > &xq, const std::string &mode, bool equidistant)
 Functions called by friend functions defined here. More...
 
static casadi_int sprank (const MatType &x)
 Functions called by friend functions defined here. More...
 
static casadi_int norm_0_mul (const MatType &x, const MatType &y)
 Functions called by friend functions defined here. More...
 
static MatType tril (const MatType &x, bool includeDiagonal=true)
 Functions called by friend functions defined here. More...
 
static MatType triu (const MatType &x, bool includeDiagonal=true)
 Functions called by friend functions defined here. More...
 
static MatType sumsqr (const MatType &x)
 Functions called by friend functions defined here. More...
 
static MatType linspace (const MatType &a, const MatType &b, casadi_int nsteps)
 Functions called by friend functions defined here. More...
 
static MatType cross (const MatType &a, const MatType &b, casadi_int dim=-1)
 Functions called by friend functions defined here. More...
 
static MatType skew (const MatType &a)
 Functions called by friend functions defined here. More...
 
static MatType inv_skew (const MatType &a)
 Functions called by friend functions defined here. More...
 
static MatType tril2symm (const MatType &x)
 Functions called by friend functions defined here. More...
 
static MatType triu2symm (const MatType &x)
 Functions called by friend functions defined here. More...
 
static MatType repsum (const MatType &x, casadi_int n, casadi_int m=1)
 Functions called by friend functions defined here. More...
 
static MatType diff (const MatType &x, casadi_int n=1, casadi_int axis=-1)
 Functions called by friend functions defined here. More...
 
static bool is_linear (const MatType &expr, const MatType &var)
 Functions called by friend functions defined here. More...
 
static bool is_quadratic (const MatType &expr, const MatType &var)
 Functions called by friend functions defined here. More...
 
static void quadratic_coeff (const MatType &expr, const MatType &var, MatType &A, MatType &b, MatType &c, bool check)
 Functions called by friend functions defined here. More...
 
static void linear_coeff (const MatType &expr, const MatType &var, MatType &A, MatType &b, bool check)
 Functions called by friend functions defined here. More...
 
static MatType jtimes (const MatType &ex, const MatType &arg, const MatType &v, bool tr=false, const Dict &opts=Dict())
 
static MatType gradient (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 
static MatType tangent (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 
static MatType linearize (const MatType &f, const MatType &x, const MatType &x0, const Dict &opts=Dict())
 
static MatType mpower (const MatType &x, const MatType &y)
 
static MatType soc (const MatType &x, const MatType &y)
 
Construct symbolic primitives

The "sym" function is intended to work in a similar way as "sym" used in the Symbolic Toolbox for Matlab but instead creating a CasADi symbolic primitive.

static MatType sym (const std::string &name, casadi_int nrow=1, casadi_int ncol=1)
 Create an nrow-by-ncol symbolic primitive. More...
 
static MatType sym (const std::string &name, const std::pair< casadi_int, casadi_int > &rc)
 Construct a symbolic primitive with given dimensions. More...
 
static MatType sym (const std::string &name, const Sparsity &sp)
 Create symbolic primitive with a given sparsity pattern. More...
 
static std::vector< MatType > sym (const std::string &name, const Sparsity &sp, casadi_int p)
 Create a vector of length p with with matrices. More...
 
static std::vector< MatType > sym (const std::string &name, casadi_int nrow, casadi_int ncol, casadi_int p)
 Create a vector of length p with nrow-by-ncol symbolic primitives. More...
 
static std::vector< std::vector< MatType > > sym (const std::string &name, const Sparsity &sp, casadi_int p, casadi_int r)
 Create a vector of length r of vectors of length p with. More...
 
static std::vector< std::vector< MatType > > sym (const std::string &name, casadi_int nrow, casadi_int ncol, casadi_int p, casadi_int r)
 Create a vector of length r of vectors of length p. More...
 
static MatType zeros (casadi_int nrow=1, casadi_int ncol=1)
 Create a dense matrix or a matrix with specified sparsity with all entries zero. More...
 
static MatType zeros (const Sparsity &sp)
 Create a dense matrix or a matrix with specified sparsity with all entries zero. More...
 
static MatType zeros (const std::pair< casadi_int, casadi_int > &rc)
 Create a dense matrix or a matrix with specified sparsity with all entries zero. More...
 
static MatType ones (casadi_int nrow=1, casadi_int ncol=1)
 Create a dense matrix or a matrix with specified sparsity with all entries one. More...
 
static MatType ones (const Sparsity &sp)
 Create a dense matrix or a matrix with specified sparsity with all entries one. More...
 
static MatType ones (const std::pair< casadi_int, casadi_int > &rc)
 Create a dense matrix or a matrix with specified sparsity with all entries one. More...
 

Friends

MatType interp1d (const std::vector< double > &x, const MatType &v, const std::vector< double > &xq, const std::string &mode, bool equidistant=false)
 Performs 1d linear interpolation. More...
 
MatType mpower (const MatType &x, const MatType &n)
 Matrix power x^n. More...
 
MatType soc (const MatType &x, const MatType &y)
 Construct second-order-convex. More...
 
MatType mrdivide (const MatType &x, const MatType &n)
 Matrix divide (cf. slash '/' in MATLAB) More...
 
MatType mldivide (const MatType &x, const MatType &n)
 Matrix divide (cf. backslash '\' in MATLAB) More...
 
std::vector< MatType > symvar (const MatType &x)
 Get all symbols contained in the supplied expression. More...
 
MatType sumsqr (const MatType &x)
 Calculate sum of squares: sum_ij X_ij^2. More...
 
MatType logsumexp (const MatType &x)
 x -> log(sum_i exp(x_i)) More...
 
MatType logsumexp (const MatType &x, const MatType &margin)
 Scaled version of logsumexp. More...
 
MatType linspace (const MatType &a, const MatType &b, casadi_int nsteps)
 Matlab's linspace command. More...
 
MatType cross (const MatType &a, const MatType &b, casadi_int dim=-1)
 Matlab's cross command. More...
 
MatType skew (const MatType &a)
 Generate a skew symmetric matrix from a 3-vector. More...
 
MatType inv_skew (const MatType &a)
 Generate the 3-vector progenitor of a skew symmetric matrix. More...
 
MatType det (const MatType &A)
 Matrix determinant (experimental) More...
 
MatType inv_minor (const MatType &A)
 Matrix inverse (experimental) More...
 
MatType inv (const MatType &A)
 Matrix inverse. More...
 
MatType inv (const MatType &A, const std::string &lsolver, const Dict &options=Dict())
 Matrix inverse. More...
 
MatType trace (const MatType &x)
 Matrix trace. More...
 
MatType tril2symm (const MatType &a)
 Convert a lower triangular matrix to a symmetric one. More...
 
MatType triu2symm (const MatType &a)
 Convert a upper triangular matrix to a symmetric one. More...
 
MatType norm_fro (const MatType &x)
 Frobenius norm. More...
 
MatType norm_2 (const MatType &x)
 2-norm More...
 
MatType norm_1 (const MatType &x)
 1-norm More...
 
MatType norm_inf (const MatType &x)
 Infinity-norm. More...
 
MatType diff (const MatType &x, casadi_int n=1, casadi_int axis=-1)
 Returns difference (n-th order) along given axis (MATLAB convention) More...
 
MatType cumsum (const MatType &x, casadi_int axis=-1)
 Returns cumulative sum along given axis (MATLAB convention) More...
 
MatType dot (const MatType &x, const MatType &y)
 Inner product of two matrices. More...
 
MatType nullspace (const MatType &A)
 Computes the nullspace of a matrix A. More...
 
MatType polyval (const MatType &p, const MatType &x)
 Evaluate a polynomial with coefficients p in x. More...
 
MatType diag (const MatType &A)
 Get the diagonal of a matrix or construct a diagonal. More...
 
MatType unite (const MatType &A, const MatType &B)
 Unite two matrices no overlapping sparsity. More...
 
MatType densify (const MatType &x)
 Make the matrix dense if not already. More...
 
MatType densify (const MatType &x, const MatType &val)
 Make the matrix dense and assign nonzeros to a value. More...
 
MatType project (const MatType &A, const Sparsity &sp, bool intersect=false)
 Create a new matrix with a given sparsity pattern but with the. More...
 
MatType if_else (const MatType &cond, const MatType &if_true, const MatType &if_false, bool short_circuit=false)
 Branching on MX nodes. More...
 
MatType conditional (const MatType &ind, const std::vector< MatType > &x, const MatType &x_default, bool short_circuit=false)
 Create a switch. More...
 
bool depends_on (const MatType &f, const MatType &arg)
 Check if expression depends on the argument. More...
 
MatType substitute (const MatType &ex, const MatType &v, const MatType &vdef)
 Substitute variable v with expression vdef in an expression ex. More...
 
std::vector< MatType > substitute (const std::vector< MatType > &ex, const std::vector< MatType > &v, const std::vector< MatType > &vdef)
 Substitute variable var with expression expr in multiple expressions. More...
 
void substitute_inplace (const std::vector< MatType > &v, std::vector< MatType > &inout_vdef, std::vector< MatType > &inout_ex, bool reverse=false)
 Inplace substitution with piggyback expressions. More...
 
MatType cse (const MatType &e)
 Common subexpression elimination. More...
 
std::vector< MatType > cse (const std::vector< MatType > &e)
 Common subexpression elimination. More...
 
MatType solve (const MatType &A, const MatType &b)
 Solve a system of equations: A*x = b. More...
 
MatType solve (const MatType &A, const MatType &b, const std::string &lsolver, const Dict &dict=Dict())
 Solve a system of equations: A*x = b. More...
 
MatType linearize (const MatType &f, const MatType &x, const MatType &x0, const Dict &opts=Dict())
 Linearize an expression. More...
 
MatType pinv (const MatType &A)
 Computes the Moore-Penrose pseudo-inverse. More...
 
MatType pinv (const MatType &A, const std::string &lsolver, const Dict &dict=Dict())
 Computes the Moore-Penrose pseudo-inverse. More...
 
MatType expm_const (const MatType &A, const MatType &t)
 Calculate Matrix exponential. More...
 
MatType expm (const MatType &A)
 Calculate Matrix exponential. More...
 
MatType jacobian (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 Calculate Jacobian. More...
 
MatType gradient (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 Calculate the gradient of an expression. More...
 
MatType tangent (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 Calculate the tangent of an expression. More...
 
MatType jtimes (const MatType &ex, const MatType &arg, const MatType &v, bool tr=false, const Dict &opts=Dict())
 Calculate the Jacobian and multiply by a vector from the right. More...
 
std::vector< std::vector< MatType > > forward (const std::vector< MatType > &ex, const std::vector< MatType > &arg, const std::vector< std::vector< MatType > > &v, const Dict &opts=Dict())
 Forward directional derivative. More...
 
std::vector< std::vector< MatType > > reverse (const std::vector< MatType > &ex, const std::vector< MatType > &arg, const std::vector< std::vector< MatType > > &v, const Dict &opts=Dict())
 Reverse directional derivative. More...
 
std::vector< bool > which_depends (const MatType &expr, const MatType &var, casadi_int order, bool tr)
 Find out which variables enter with some order. More...
 
Sparsity jacobian_sparsity (const MatType &f, const MatType &x)
 Get the sparsity pattern of a jacobian. More...
 
bool is_linear (const MatType &expr, const MatType &var)
 Is expr linear in var? More...
 
bool is_quadratic (const MatType &expr, const MatType &var)
 Is expr quadratic in var? More...
 
void quadratic_coeff (const MatType &expr, const MatType &var, MatType &A, MatType &b, MatType &c, bool check=true)
 Recognizes quadratic form in scalar expression. More...
 
void linear_coeff (const MatType &expr, const MatType &var, MatType &A, MatType &b, bool check=true)
 Recognizes linear form in vector expression. More...
 
void extract_parametric (const MatType &expr, const MatType &par, MatType &expr_ret, std::vector< MatType > &symbols, std::vector< MatType > &parametric, const Dict &opts=Dict())
 Extract purely parametric parts from an expression graph. More...
 
void extract_parametric (const std::vector< MatType > &expr, const MatType &par, std::vector< MatType > &expr_ret, std::vector< MatType > &symbols, std::vector< MatType > &parametric, const Dict &opts=Dict())
 
void extract_parametric (const std::vector< MatType > &expr, const std::vector< MatType > &par, std::vector< MatType > &expr_ret, std::vector< MatType > &symbols, std::vector< MatType > &parametric, const Dict &opts=Dict())
 
void extract_parametric (const MatType &expr, const std::vector< MatType > &par, MatType &expr_ret, std::vector< MatType > &symbols, std::vector< MatType > &parametric, const Dict &opts=Dict())
 
void separate_linear (const MatType &expr, const MatType &sym_lin, const MatType &sym_const, MatType &expr_const, MatType &expr_lin, MatType &expr_nonlin)
 
void separate_linear (const MatType &expr, const std::vector< MatType > &sym_lin, const std::vector< MatType > &sym_const, MatType &expr_const, MatType &expr_lin, MatType &expr_nonlin)
 
casadi_int n_nodes (const MatType &A)
 
MatType simplify (const MatType &x)
 Simplify an expression. More...
 
std::string print_operator (const MatType &xb, const std::vector< std::string > &args)
 Get a string representation for a binary MatType, using custom arguments. More...
 
void extract (std::vector< MatType > &ex, std::vector< MatType > &v, std::vector< MatType > &vdef, const Dict &opts=Dict())
 Introduce intermediate variables for selected nodes in a graph. More...
 
void shared (std::vector< MatType > &ex, std::vector< MatType > &v, std::vector< MatType > &vdef, const std::string &v_prefix="v_", const std::string &v_suffix="")
 Extract shared subexpressions from an set of expressions. More...
 
MatType repsum (const MatType &A, casadi_int n, casadi_int m=1)
 Given a repeated matrix, computes the sum of repeated parts. More...
 
MatType einstein (const MatType &A, const MatType &B, const MatType &C, const std::vector< casadi_int > &dim_a, const std::vector< casadi_int > &dim_b, const std::vector< casadi_int > &dim_c, const std::vector< casadi_int > &a, const std::vector< casadi_int > &b, const std::vector< casadi_int > &c)
 Compute any contraction of two dense tensors, using index/einstein notation. More...
 
MatType einstein (const MatType &A, const MatType &B, const std::vector< casadi_int > &dim_a, const std::vector< casadi_int > &dim_b, const std::vector< casadi_int > &dim_c, const std::vector< casadi_int > &a, const std::vector< casadi_int > &b, const std::vector< casadi_int > &c)
 Compute any contraction of two dense tensors, using index/einstein notation. More...
 
bool contains (const std::vector< MatType > &v, const MatType &n)
 Check if expression n is listed in v. More...
 
bool contains_all (const std::vector< MatType > &v, const std::vector< MatType > &n)
 Check if expression n is listed in v. More...
 
bool contains_any (const std::vector< MatType > &v, const std::vector< MatType > &n)
 Check if expression n is listed in v. More...
 
MatType hessian (const MatType &ex, const MatType &arg, const Dict &opts=Dict())
 Hessian and (optionally) gradient. More...
 
MatType hessian (const MatType &ex, const MatType &arg, MatType &output_g, const Dict &opts=Dict())
 Hessian and (optionally) gradient. More...
 
MatType mmin (const MatType &x)
 Smallest element in a matrix. More...
 
MatType mmax (const MatType &x)
 Largest element in a matrix. More...
 
static MatType bilin (const MatType &A, const MatType &x, const MatType &y)
 Calculate bilinear/quadratic form x^T A y. More...
 
MatType bilin (const MatType &A, const MatType &x, const MatType &y)
 Calculate bilinear/quadratic form x^T A y. More...
 
MatType bilin (const MatType &A, const MatType &x)
 Calculate bilinear/quadratic form x^T A y. More...
 
static MatType rank1 (const MatType &A, const MatType &alpha, const MatType &x, const MatType &y)
 Make a rank-1 update to a matrix A. More...
 
MatType rank1 (const MatType &A, const MatType &alpha, const MatType &x, const MatType &y)
 Make a rank-1 update to a matrix A. More...
 

Member Function Documentation

◆ colind() [1/2]

template<typename MatType >
const casadi_int* casadi::GenericMatrix< MatType >::colind ( ) const
inline

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

Definition at line 198 of file generic_matrix.hpp.

198 { return sparsity().colind(); }
Sparsity sparsity() const
Get the sparsity pattern.
const casadi_int * colind() const
Get a reference to the colindex of all column element (see class description)
Definition: sparsity.cpp:168

References casadi::Sparsity::colind(), and casadi::GenericMatrix< MatType >::sparsity().

Referenced by casadi::Matrix< Scalar >::densify(), casadi::SnoptInterface::init(), casadi::Conic::sdp_to_socp_init(), casadi::MX::set(), casadi::SnoptInterface::solve(), and casadi::SnoptInterface::userfun().

◆ colind() [2/2]

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::colind ( casadi_int  col) const
inline

◆ columns()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::columns ( ) const
inline

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

Definition at line 124 of file generic_matrix.hpp.

124 {return size2();}
casadi_int size2() const
Get the second dimension (i.e. number of columns)

References casadi::GenericMatrix< MatType >::size2().

◆ cross()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::cross ( const MatType &  a,
const MatType &  b,
casadi_int  dim = -1 
)
static

Definition at line 1394 of file generic_matrix.hpp.

1394  {
1395  casadi_assert(a.size1()==b.size1() && a.size2()==b.size2(),
1396  "cross(a, b): Inconsistent dimensions. Dimension of a ("
1397  + a.dim() + " ) must equal that of b (" + b.dim() + ").");
1398 
1399  casadi_assert(a.size1()==3 || a.size2()==3,
1400  "cross(a, b): One of the dimensions of a should have length 3, but got "
1401  + a.dim() + ".");
1402  casadi_assert(dim==-1 || dim==1 || dim==2,
1403  "cross(a, b, dim): Dim must be 1, 2 or -1 (automatic).");
1404 
1405  std::vector<MatType> ret(3);
1406 
1407  bool t = a.size1()==3;
1408 
1409  if (dim==1) t = true;
1410  if (dim==2) t = false;
1411 
1412  MatType a1 = t ? a(0, Slice()) : a(Slice(), 0);
1413  MatType a2 = t ? a(1, Slice()) : a(Slice(), 1);
1414  MatType a3 = t ? a(2, Slice()) : a(Slice(), 2);
1415 
1416  MatType b1 = t ? b(0, Slice()) : b(Slice(), 0);
1417  MatType b2 = t ? b(1, Slice()) : b(Slice(), 1);
1418  MatType b3 = t ? b(2, Slice()) : b(Slice(), 2);
1419 
1420  ret[0] = a2*b3-a3*b2;
1421  ret[1] = a3*b1-a1*b3;
1422  ret[2] = a1*b2-a2*b1;
1423 
1424  return t ? vertcat(ret) : horzcat(ret);
1425  }
std::string dim(bool with_nz=false) const
Get string representation of dimensions.
friend MatType vertcat(const std::vector< MatType > &v)
Concatenate a list of matrices vertically.
friend MatType horzcat(const std::vector< MatType > &v)
Concatenate a list of matrices horizontally.

◆ diff()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::diff ( const MatType &  x,
casadi_int  n = 1,
casadi_int  axis = -1 
)
static

Definition at line 1774 of file generic_matrix.hpp.

1774  {
1775  casadi_assert(axis==-1 || axis==0 || axis==1, "Axis argument invalid");
1776  casadi_assert(n>=1, "n argument invalid");
1777 
1778  MatType ret = x;
1779  for (casadi_int i=0;i<n;++i) {
1780  // Matlab's special case
1781  if (axis==-1 && ret.is_scalar()) return MatType();
1782 
1783  casadi_int local_axis = (axis==-1) ? ret.is_row() : axis;
1784  if (local_axis==0) {
1785  if (ret.size1()<=1) {
1786  ret = MatType::zeros(0, ret.size2());
1787  } else {
1788  ret = ret(Slice(1, ret.size1()), Slice())-ret(Slice(0, ret.size1()-1), Slice());
1789  }
1790  } else {
1791  if (ret.size2()<=1) {
1792  ret = MatType::zeros(ret.size1(), 0);
1793  } else {
1794  ret = ret(Slice(), Slice(1, ret.size2()))-ret(Slice(), Slice(0, ret.size2()-1));
1795  }
1796  }
1797  }
1798  return ret;
1799  }

◆ dim()

template<typename MatType >
std::string casadi::GenericMatrix< MatType >::dim ( bool  with_nz = false) const

◆ get_colind()

template<typename MatType >
std::vector<casadi_int> casadi::GenericMatrix< MatType >::get_colind ( ) const
inline

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

Definition at line 195 of file generic_matrix.hpp.

195 { return sparsity().get_colind(); }
std::vector< casadi_int > get_colind() const
Get the column index for each column.
Definition: sparsity.cpp:364

References casadi::Sparsity::get_colind(), and casadi::GenericMatrix< MatType >::sparsity().

◆ get_row()

template<typename MatType >
std::vector<casadi_int> casadi::GenericMatrix< MatType >::get_row ( ) const
inline

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

Definition at line 194 of file generic_matrix.hpp.

194 { return sparsity().get_row(); }
std::vector< casadi_int > get_row() const
Get the row for each non-zero entry.
Definition: sparsity.cpp:372

References casadi::Sparsity::get_row(), and casadi::GenericMatrix< MatType >::sparsity().

◆ interp1d()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::interp1d ( const std::vector< double > &  x,
const MatType &  v,
const std::vector< double > &  xq,
const std::string &  mode,
bool  equidistant 
)
static

Definition at line 1431 of file generic_matrix.hpp.

1432  {
1433 
1434  bool mode_floor = false;
1435  bool mode_ceil = false;
1436  if (mode=="floor") {
1437  mode_floor = true;
1438  } else if (mode=="ceil") {
1439  mode_ceil = true;
1440  } else if (mode=="linear") {
1441  //
1442  } else {
1443  casadi_error("interp1d(x, v, xq, mode): "
1444  "Mode must be 'floor', 'ceil' or 'linear'. Got '" + mode + "' instead.");
1445  }
1446 
1447  casadi_assert(is_increasing(x), "interp1d(x, v, xq): x must be increasing.");
1448 
1449  casadi_assert(x.size()==v.size1(),
1450  "interp1d(x, v, xq): dimensions mismatch. v expected to have " + str(x.size()) + " rows,"
1451  " but got " + str(v.size1()) + " instead.");
1452 
1453  // Need at least two elements
1454  casadi_assert(x.size()>=2, "interp1d(x, v, xq): x must be at least length 2.");
1455 
1456  // Vectors to compose a sparse matrix
1457  std::vector<double> val;
1458  std::vector<casadi_int> colind(1, 0);
1459  std::vector<casadi_int> row;
1460 
1461  // Number of nonzeros in to-be composed matrix
1462  casadi_int nnz = 0;
1463  for (casadi_int i=0;i<xq.size();++i) {
1464  // Obtain index corresponding to xq[i]
1465  double ind = index_interp1d(x, xq[i], equidistant);
1466 
1467  if (mode_floor) ind = floor(ind);
1468  if (mode_ceil) ind = ceil(ind);
1469 
1470  // Split into integer and fractional part
1471  double int_partd;
1472  double frac_part = modf(ind, &int_partd);
1473  casadi_int int_part = static_cast<casadi_int>(int_partd);
1474 
1475  if (frac_part==0) {
1476  // Create a single entry
1477  val.push_back(1);
1478  row.push_back(int_part);
1479  nnz+=1;
1480  colind.push_back(nnz);
1481  } else {
1482  // Create a double entry
1483  val.push_back(1-frac_part);
1484  val.push_back(frac_part);
1485  row.push_back(int_part);
1486  row.push_back(int_part+1);
1487  nnz+=2;
1488  colind.push_back(nnz);
1489  }
1490  }
1491 
1492  // Construct sparsity for composed matrix
1493  Sparsity sp(x.size(), xq.size() , colind, row);
1494 
1495  return MatType::mtimes(MatType(sp, val).T(), v);
1496 
1497  }
casadi_int nnz() const
Get the number of (structural) non-zero elements.
const casadi_int * colind() const
Get the sparsity pattern. See the Sparsity class for details.
const casadi_int * row() const
Get the sparsity pattern. See the Sparsity class for details.
double index_interp1d(const std::vector< double > &x, double xq, bool equidistant)
bool is_increasing(const std::vector< T > &v)
Check if the vector is strictly increasing.
std::string str(const T &v)
String representation, any type.

References casadi::index_interp1d(), casadi::is_increasing(), casadi::str(), and casadi::T.

◆ inv_skew()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::inv_skew ( const MatType &  a)
static

Definition at line 1511 of file generic_matrix.hpp.

1511  {
1512  casadi_assert(a.size1()==3 && a.size2()==3,
1513  "inv_skew(a): Expecting 3-by-3 matrix, got " + a.dim() + ".");
1514 
1515  return 0.5*vertcat(std::vector<MatType>({a(2, 1)-a(1, 2), a(0, 2)-a(2, 0), a(1, 0)-a(0, 1)}));
1516  }

◆ is_column()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_column ( ) const
inline

◆ is_dense()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_dense ( ) const
inline

◆ is_empty()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_empty ( bool  both = false) const
inline

◆ is_linear()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_linear ( const MatType &  expr,
const MatType &  var 
)
static

Definition at line 1743 of file generic_matrix.hpp.

1743  {
1744  return !any(MatType::which_depends(expr, var, 2, true));
1745  }
bool any(const std::vector< bool > &v)
Check if any arguments are true.
Definition: casadi_misc.cpp:84

References casadi::any().

◆ is_quadratic()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_quadratic ( const MatType &  expr,
const MatType &  var 
)
static

Definition at line 1748 of file generic_matrix.hpp.

1748  {
1749  return is_linear(gradient(expr, var), var);
1750  }
static bool is_linear(const MatType &expr, const MatType &var)
Functions called by friend functions defined here.
friend MatType gradient(const MatType &ex, const MatType &arg, const Dict &opts=Dict())
Calculate the gradient of an expression.

◆ is_row()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_row ( ) const
inline

◆ is_scalar()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_scalar ( bool  scalar_and_dense = false) const

◆ is_square()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_square ( ) const
inline

◆ is_tril()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_tril ( ) const
inline

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

Definition at line 188 of file generic_matrix.hpp.

188 { return sparsity().is_tril();}
bool is_tril(bool strictly=false) const
Is lower triangular?
Definition: sparsity.cpp:321

References casadi::Sparsity::is_tril(), and casadi::GenericMatrix< MatType >::sparsity().

Referenced by casadi::MX::solve().

◆ is_triu()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_triu ( ) const
inline

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

Definition at line 183 of file generic_matrix.hpp.

183 { return sparsity().is_triu();}
bool is_triu(bool strictly=false) const
Is upper triangular?
Definition: sparsity.cpp:325

References casadi::Sparsity::is_triu(), and casadi::GenericMatrix< MatType >::sparsity().

Referenced by casadi::MX::solve().

◆ is_vector()

template<typename MatType >
bool casadi::GenericMatrix< MatType >::is_vector ( ) const
inline

◆ linear_coeff()

template<typename MatType >
void casadi::GenericMatrix< MatType >::linear_coeff ( const MatType &  expr,
const MatType &  var,
MatType &  A,
MatType &  b,
bool  check 
)
static

Definition at line 1764 of file generic_matrix.hpp.

1765  {
1766  casadi_assert(expr.is_vector(), "'linear_coeff' only defined for vector expressions.");
1767  if (check)
1768  casadi_assert(is_linear(expr, var), "'linear_coeff' called on non-linear expression.");
1769  A = substitute(jacobian(expr, var), var, 0);
1770  b = vec(substitute(expr, var, 0));
1771  }
static MatType vec(const MatType &x)
friend MatType substitute(const MatType &ex, const MatType &v, const MatType &vdef)
Substitute variable v with expression vdef in an expression ex.
friend MatType jacobian(const MatType &ex, const MatType &arg, const Dict &opts=Dict())
Calculate Jacobian.

◆ linspace()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::linspace ( const MatType &  a,
const MatType &  b,
casadi_int  nsteps 
)
static

Definition at line 1381 of file generic_matrix.hpp.

1381  {
1382  std::vector<MatType> ret(nsteps);
1383  ret[0] = a;
1384  MatType step = (b-a)/static_cast<MatType>(nsteps-1);
1385 
1386  for (casadi_int i=1; i<nsteps-1; ++i)
1387  ret[i] = a + i * step;
1388 
1389  ret[nsteps-1] = b;
1390  return vertcat(ret);
1391  }

◆ nnz()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::nnz

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

Definition at line 1298 of file generic_matrix.hpp.

1298  {
1299  return sparsity().nnz();
1300  }
casadi_int nnz() const
Get the number of (structural) non-zeros.
Definition: sparsity.cpp:148

Referenced by casadi::Constant< Value >::_get_binary(), casadi::OptiNode::bake(), casadi::ConstantMX::create(), casadi::GetNonzerosParam::create(), casadi::SetNonzeros< Add >::create(), casadi::GetNonzeros::create(), casadi::MX::depends_on(), casadi::Diagsplit::Diagsplit(), casadi::Find::eval(), casadi::Low::eval(), casadi::GetNonzerosParamVector::eval(), casadi::GetNonzerosSliceParam::eval(), casadi::GetNonzerosParamSlice::eval(), casadi::GetNonzerosParamParam::eval(), casadi::HorzRepmat::eval_gen(), casadi::Concat::eval_gen(), casadi::Transpose::eval_gen(), casadi::Find::generate(), casadi::Low::generate(), casadi::Concat::generate(), casadi::GetNonzerosParamVector::generate(), casadi::GetNonzerosSliceParam::generate(), casadi::GetNonzerosParamSlice::generate(), casadi::GetNonzerosParamParam::generate(), casadi::Output::generate(), casadi::HorzRepmat::generate(), casadi::MXNode::get_binary(), casadi::MX::get_nonzeros(), casadi::MX::get_nz(), casadi::Concat::get_nzref(), casadi::AmplInterface::init(), casadi::SnoptInterface::init(), casadi::SuperscsInterface::init(), casadi::Scpgen::init(), casadi::CplexInterface::init_mem(), casadi::SnoptInterface::init_mem(), casadi::SymbolicQr::linsol_eval_sx(), casadi::Monitor::Monitor(), casadi::mtaylor_recursive(), casadi::MX::polyval(), casadi::Reshape::Reshape(), casadi::Matrix< Scalar >::set(), casadi::MX::set(), casadi::MX::set_nz(), casadi::SuperscsInterface::solve(), casadi::MX::solve(), casadi::SnoptInterface::solve(), casadi::MXNode::sp_forward(), casadi::Concat::sp_forward(), casadi::Dot::sp_forward(), casadi::MXNode::sp_reverse(), casadi::Concat::sp_reverse(), casadi::Dot::sp_reverse(), casadi::GetNonzerosParam::sp_reverse(), casadi::HorzRepmat::sp_reverse(), casadi::MX::sparsity_cast(), casadi::SparsityCast::SparsityCast(), casadi::OptiNode::subject_to(), casadi::GetNonzerosParamSlice::sz_iw(), and casadi::GetNonzerosParamParam::sz_iw().

◆ nnz_diag()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::nnz_diag

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

Definition at line 1313 of file generic_matrix.hpp.

1313  {
1314  return sparsity().nnz_diag();
1315  }
casadi_int nnz_diag() const
Number of non-zeros on the diagonal, i.e. the number of elements (i, j) with j==i.
Definition: sparsity.cpp:360

◆ nnz_lower()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::nnz_lower

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

Definition at line 1303 of file generic_matrix.hpp.

1303  {
1304  return sparsity().nnz_lower();
1305  }
casadi_int nnz_lower(bool strictly=false) const
Number of non-zeros in the lower triangular half,.
Definition: sparsity.cpp:352

◆ nnz_upper()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::nnz_upper

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

Definition at line 1308 of file generic_matrix.hpp.

1308  {
1309  return sparsity().nnz_upper();
1310  }
casadi_int nnz_upper(bool strictly=false) const
Number of non-zeros in the upper triangular half,.
Definition: sparsity.cpp:356

◆ norm_0_mul()

template<typename MatType >
static casadi_int casadi::GenericMatrix< MatType >::norm_0_mul ( const MatType &  x,
const MatType &  y 
)
inlinestatic

Definition at line 216 of file generic_matrix.hpp.

216  {
217  return Sparsity::norm_0_mul(x.sparsity(), y.sparsity());
218  }
static casadi_int norm_0_mul(const Sparsity &x, const Sparsity &A)
Enlarge matrix.
Definition: sparsity.cpp:1665

References casadi::Sparsity::norm_0_mul().

◆ numel()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::numel

◆ nz() [1/2]

template<typename MatType >
template<typename K >
NonZeros<MatType, K> casadi::GenericMatrix< MatType >::nz ( const K &  k)
inline

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

Definition at line 258 of file generic_matrix.hpp.

258  {
259  return NonZeros<MatType, K>(self(), k);
260  }

◆ nz() [2/2]

template<typename MatType >
template<typename K >
const MatType casadi::GenericMatrix< MatType >::nz ( const K &  k) const
inline

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

Definition at line 248 of file generic_matrix.hpp.

248  {
249  MatType ret;
250  self().get_nz(ret, false, k);
251  return ret;
252  }

Referenced by casadi::SnoptInterface::init(), casadi::MX::kron(), and casadi::MX::polyval().

◆ ones() [1/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::ones ( casadi_int  nrow = 1,
casadi_int  ncol = 1 
)
inlinestatic

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

Definition at line 1275 of file generic_matrix.hpp.

1275  {
1276  return ones(Sparsity::dense(nrow, ncol));
1277  }
static MatType ones(casadi_int nrow=1, casadi_int ncol=1)
Create a dense matrix or a matrix with specified sparsity with all entries one.
static Sparsity dense(casadi_int nrow, casadi_int ncol=1)
Create a dense rectangular sparsity pattern *.
Definition: sparsity.cpp:1012

References casadi::Sparsity::dense().

Referenced by casadi::OptiNode::bake(), casadi::OptiNode::canon_expr(), casadi::GenericMatrix< MatType >::ones(), and casadi::OptiNode::variable().

◆ ones() [2/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::ones ( const Sparsity sp)
inlinestatic

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

Definition at line 1278 of file generic_matrix.hpp.

1278 { return MatType(sp, 1, false);}

◆ ones() [3/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::ones ( const std::pair< casadi_int, casadi_int > &  rc)
inlinestatic

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

Definition at line 1279 of file generic_matrix.hpp.

1279  {
1280  return ones(rc.first, rc.second);
1281  }

References casadi::GenericMatrix< MatType >::ones().

◆ operator()() [1/4]

template<typename MatType >
template<typename RR >
SubIndex<MatType, RR> casadi::GenericMatrix< MatType >::operator() ( const RR &  rr)
inline

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

Definition at line 286 of file generic_matrix.hpp.

286  {
287  return SubIndex<MatType, RR>(self(), rr);
288  }

◆ operator()() [2/4]

template<typename MatType >
template<typename RR >
const MatType casadi::GenericMatrix< MatType >::operator() ( const RR &  rr) const
inline

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

Definition at line 266 of file generic_matrix.hpp.

266  {
267  MatType ret;
268  self().get(ret, false, rr);
269  return ret;
270  }

◆ operator()() [3/4]

template<typename MatType >
template<typename RR , typename CC >
SubMatrix<MatType, RR, CC> casadi::GenericMatrix< MatType >::operator() ( const RR &  rr,
const CC &  cc 
)
inline

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

Definition at line 294 of file generic_matrix.hpp.

294  {
295  return SubMatrix<MatType, RR, CC>(self(), rr, cc);
296  }

◆ operator()() [4/4]

template<typename MatType >
template<typename RR , typename CC >
const MatType casadi::GenericMatrix< MatType >::operator() ( const RR &  rr,
const CC &  cc 
) const
inline

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

Definition at line 276 of file generic_matrix.hpp.

276  {
277  MatType ret;
278  self().get(ret, false, rr, cc);
279  return ret;
280  }

◆ quadratic_coeff()

template<typename MatType >
void casadi::GenericMatrix< MatType >::quadratic_coeff ( const MatType &  expr,
const MatType &  var,
MatType &  A,
MatType &  b,
MatType &  c,
bool  check 
)
static

Definition at line 1753 of file generic_matrix.hpp.

1754  {
1755  casadi_assert(expr.is_scalar(), "'quadratic_coeff' only defined for scalar expressions.");
1756  A = hessian(expr, var);
1757  b = substitute(jacobian(expr, var), var, 0).T();
1758  if (check)
1759  casadi_assert(!depends_on(A, var), "'quadratic_coeff' called on non-quadratic expression.");
1760  c = substitute(expr, var, 0);
1761  }
friend bool depends_on(const MatType &f, const MatType &arg)
Check if expression depends on the argument.
friend MatType hessian(const MatType &ex, const MatType &arg, const Dict &opts=Dict())
Hessian and (optionally) gradient.

◆ repsum()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::repsum ( const MatType &  x,
casadi_int  n,
casadi_int  m = 1 
)
static

Definition at line 1529 of file generic_matrix.hpp.

1529  {
1530  casadi_assert_dev(x.size1() % n==0);
1531  casadi_assert_dev(x.size2() % m==0);
1532  std::vector< std::vector< MatType> > s =
1533  blocksplit(x, x.size1()/n, x.size2()/m);
1534  MatType sum = 0;
1535  for (casadi_int i=0;i<s.size();++i) {
1536  for (casadi_int j=0;j<s[i].size();++j) {
1537  sum = sum + s[i][j];
1538  }
1539  }
1540  return sum;
1541  }
static std::vector< std::vector< MatType > > blocksplit(const MatType &x, const std::vector< casadi_int > &vert_offset, const std::vector< casadi_int > &horz_offset)
friend MatType sum(const MatType &x)
Returns summation of all elements.

References casadi::sum().

Referenced by casadi::MXNode::get_repsum().

◆ row() [1/2]

template<typename MatType >
const casadi_int* casadi::GenericMatrix< MatType >::row ( ) const
inline

◆ row() [2/2]

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::row ( casadi_int  el) const
inline

◆ rows()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::rows ( ) const
inline

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

Definition at line 114 of file generic_matrix.hpp.

114 {return size1();}
casadi_int size1() const
Get the first dimension (i.e. number of rows)

References casadi::GenericMatrix< MatType >::size1().

◆ size() [1/2]

template<typename MatType >
std::pair< casadi_int, casadi_int > casadi::GenericMatrix< MatType >::size

◆ size() [2/2]

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::size ( casadi_int  axis) const

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

Definition at line 1338 of file generic_matrix.hpp.

1338  {
1339  return sparsity().size(axis);
1340  }

◆ size1()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::size1

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

Definition at line 1323 of file generic_matrix.hpp.

1323  {
1324  return sparsity().size1();
1325  }
casadi_int size1() const
Get the number of rows.
Definition: sparsity.cpp:124

Referenced by casadi::Constant< Value >::_get_binary(), casadi::MX::binary(), casadi::BSplineParametric::BSplineParametric(), casadi::OptiNode::canon_expr(), casadi::BSplineInterpolant::construct_graph(), casadi::Matrix< Scalar >::cumsum(), casadi::Matrix< Scalar >::densify(), casadi::MX::diagsplit(), casadi::do_inline(), casadi::dplesol(), casadi::DaeBuilderInternal::eq(), casadi::Find::eval(), casadi::DenseTranspose::eval_gen(), casadi::Find::generate(), casadi::DenseMultiplication::generate(), casadi::DenseTranspose::generate(), casadi::get_boor(), casadi::MXNode::get_dot(), casadi::MXNode::get_mac(), casadi::MXNode::get_vertsplit(), casadi::SuperscsInterface::init(), casadi::SuperscsInterface::init_mem(), casadi::MX::inv(), casadi::Inverse::Inverse(), casadi::MX::kron(), casadi::SymbolicQr::linsol_eval_sx(), casadi::MX::mac(), casadi::Multiplication::Multiplication(), casadi::MX::MX(), casadi::OptiNode::parameter(), casadi::MX::pinv(), casadi::Call::projectArg(), casadi::MX::repmat(), casadi::MX::reshape(), casadi::GenericMatrix< MatType >::rows(), casadi::Matrix< Scalar >::set(), casadi::MX::set(), casadi::Matrix< Scalar >::set_nz(), casadi::MX::set_nz(), casadi::simpleIRK(), casadi::Linsol::solve(), casadi::SuperscsInterface::solve(), casadi::Solve< Tr >::Solve(), casadi::DenseTranspose::sp_forward(), casadi::DenseTranspose::sp_reverse(), casadi::MX::sum1(), casadi::OptiNode::variable(), and casadi::MX::vertsplit().

◆ size2()

template<typename MatType >
casadi_int casadi::GenericMatrix< MatType >::size2

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

Definition at line 1328 of file generic_matrix.hpp.

1328  {
1329  return sparsity().size2();
1330  }
casadi_int size2() const
Get the number of columns.
Definition: sparsity.cpp:128

Referenced by casadi::Constant< Value >::_get_binary(), casadi::MX::binary(), casadi::OptiNode::canon_expr(), casadi::GenericMatrix< MatType >::columns(), casadi::BSplineInterpolant::construct_graph(), casadi::Matrix< Scalar >::cumsum(), casadi::MX::cumsum(), casadi::Matrix< Scalar >::densify(), casadi::MX::diagsplit(), casadi::do_inline(), casadi::DenseTranspose::eval_gen(), casadi::DenseMultiplication::generate(), casadi::DenseTranspose::generate(), casadi::get_boor(), casadi::MXNode::get_dot(), casadi::MXNode::get_horzsplit(), casadi::MXNode::get_mac(), casadi::HorzRepsum::HorzRepsum(), casadi::MX::horzsplit(), casadi::CplexInterface::init_mem(), casadi::SnoptInterface::init_mem(), casadi::SuperscsInterface::init_mem(), casadi::Inverse::Inverse(), casadi::MX::kron(), casadi::MX::mac(), casadi::Multiplication::Multiplication(), casadi::CsparseInterface::nfact(), casadi::OptiNode::parameter(), casadi::MX::pinv(), casadi::Call::projectArg(), casadi::MX::repmat(), casadi::MX::reshape(), casadi::Conic::sdp_to_socp_init(), casadi::Matrix< Scalar >::set(), casadi::MX::set(), casadi::Matrix< Scalar >::set_nz(), casadi::MX::set_nz(), casadi::Linsol::solve(), casadi::SuperscsInterface::solve(), casadi::Solve< Tr >::Solve(), casadi::DenseTranspose::sp_forward(), casadi::DenseTranspose::sp_reverse(), casadi::MX::sum2(), casadi::MX::trace(), and casadi::OptiNode::variable().

◆ skew()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::skew ( const MatType &  a)
static

Definition at line 1500 of file generic_matrix.hpp.

1500  {
1501  casadi_assert(a.is_vector() && (a.size1()==3 || a.size2()==3),
1502  "skew(a): Expecting 3-vector, got " + a.dim() + ".");
1503 
1504  MatType x = a(0);
1505  MatType y = a(1);
1506  MatType z = a(2);
1507  return blockcat(std::vector< std::vector<MatType> >({{0, -z, y}, {z, 0, -x}, {-y, x, 0}}));
1508  }
friend MatType blockcat(const std::vector< std::vector< MatType > > &v)
Construct a matrix from a list of list of blocks.

◆ sparsity()

template<typename MatType >
const Sparsity & casadi::GenericMatrix< MatType >::sparsity

◆ sprank()

template<typename MatType >
static casadi_int casadi::GenericMatrix< MatType >::sprank ( const MatType &  x)
inlinestatic

Definition at line 215 of file generic_matrix.hpp.

215 { return Sparsity::sprank(x.sparsity());}
static casadi_int sprank(const Sparsity &x)
Enlarge matrix.
Definition: sparsity.cpp:1655

References casadi::Sparsity::sprank().

◆ sumsqr()

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::sumsqr ( const MatType &  x)
inlinestatic

Definition at line 225 of file generic_matrix.hpp.

225 { return dot(x, x);}
friend MatType dot(const MatType &x, const MatType &y)
Inner product of two matrices.

References casadi::GenericMatrix< MatType >::dot.

◆ sym() [1/7]

template<typename MatType >
static std::vector<MatType > casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
casadi_int  nrow,
casadi_int  ncol,
casadi_int  p 
)
inlinestatic

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

Definition at line 1234 of file generic_matrix.hpp.

1235  {
1236  return sym(name, Sparsity::dense(nrow, ncol), p);
1237  }
static MatType sym(const std::string &name, casadi_int nrow=1, casadi_int ncol=1)
Create an nrow-by-ncol symbolic primitive.

References casadi::Sparsity::dense(), and casadi::GenericMatrix< MatType >::sym().

◆ sym() [2/7]

template<typename MatType >
static std::vector<std::vector<MatType> > casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
casadi_int  nrow,
casadi_int  ncol,
casadi_int  p,
casadi_int  r 
)
inlinestatic

with nrow-by-ncol symbolic primitives

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

Definition at line 1253 of file generic_matrix.hpp.

1253  {
1254  return sym(name, Sparsity::dense(nrow, ncol), p, r);
1255  }

References casadi::Sparsity::dense(), and casadi::GenericMatrix< MatType >::sym().

◆ sym() [3/7]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
casadi_int  nrow = 1,
casadi_int  ncol = 1 
)
inlinestatic

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

Definition at line 1206 of file generic_matrix.hpp.

1206  {
1207  return sym(name, Sparsity::dense(nrow, ncol));
1208  }

References casadi::Sparsity::dense().

Referenced by casadi::GenericMatrix< MatType >::sym().

◆ sym() [4/7]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
const Sparsity sp 
)
inlinestatic

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

Definition at line 1220 of file generic_matrix.hpp.

1220  {
1221  return MatType::_sym(name, sp);
1222  }

◆ sym() [5/7]

template<typename MatType >
std::vector< MatType > casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
const Sparsity sp,
casadi_int  p 
)
static

with symbolic primitives of given sparsity

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

Definition at line 1355 of file generic_matrix.hpp.

1356  {
1357  std::vector<MatType> ret(p);
1358  std::stringstream ss;
1359  for (casadi_int k=0; k<p; ++k) {
1360  ss.str("");
1361  ss << name << k;
1362  ret[k] = sym(ss.str(), sp);
1363  }
1364  return ret;
1365  }

◆ sym() [6/7]

template<typename MatType >
std::vector< std::vector< MatType > > casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
const Sparsity sp,
casadi_int  p,
casadi_int  r 
)
static

symbolic primitives with given sparsity

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

Definition at line 1368 of file generic_matrix.hpp.

1370  {
1371  std::vector<std::vector<MatType> > ret(r);
1372  for (casadi_int k=0; k<r; ++k) {
1373  std::stringstream ss;
1374  ss << name << "_" << k;
1375  ret[k] = sym(ss.str(), sp, p);
1376  }
1377  return ret;
1378  }

◆ sym() [7/7]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::sym ( const std::string &  name,
const std::pair< casadi_int, casadi_int > &  rc 
)
inlinestatic

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

Definition at line 1213 of file generic_matrix.hpp.

1213  {
1214  return sym(name, rc.first, rc.second);
1215  }

References casadi::GenericMatrix< MatType >::sym().

◆ tril()

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::tril ( const MatType &  x,
bool  includeDiagonal = true 
)
inlinestatic

Definition at line 219 of file generic_matrix.hpp.

219  {
220  return project(x, Sparsity::tril(x.sparsity(), includeDiagonal));
221  }
static Sparsity tril(const Sparsity &x, bool includeDiagonal=true)
Enlarge matrix.
Definition: sparsity.cpp:979
friend MatType project(const MatType &A, const Sparsity &sp, bool intersect=false)
Create a new matrix with a given sparsity pattern but with the.

References casadi::GenericMatrix< MatType >::project, and casadi::Sparsity::tril().

◆ tril2symm()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::tril2symm ( const MatType &  x)
static

Definition at line 1520 of file generic_matrix.hpp.

1520  {
1521  casadi_assert(x.is_square(),
1522  "Shape error in tril2symm. Expecting square shape but got " + x.dim());
1523  casadi_assert(x.nnz_upper()-x.nnz_diag()==0,
1524  "Sparsity error in tril2symm. Found above-diagonal entries in argument: " + x.dim());
1525  return x + x.T() - diag(diag(x));
1526  }
friend MatType diag(const MatType &A)
Get the diagonal of a matrix or construct a diagonal.

◆ triu()

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::triu ( const MatType &  x,
bool  includeDiagonal = true 
)
inlinestatic

Definition at line 222 of file generic_matrix.hpp.

222  {
223  return project(x, Sparsity::triu(x.sparsity(), includeDiagonal));
224  }
static Sparsity triu(const Sparsity &x, bool includeDiagonal=true)
Enlarge matrix.
Definition: sparsity.cpp:983

References casadi::GenericMatrix< MatType >::project, and casadi::Sparsity::triu().

◆ triu2symm()

template<typename MatType >
MatType casadi::GenericMatrix< MatType >::triu2symm ( const MatType &  x)
static

Definition at line 1544 of file generic_matrix.hpp.

1544  {
1545  casadi_assert(x.is_square(),
1546  "Shape error in triu2symm. Expecting square shape but got " + x.dim());
1547  casadi_assert(x.nnz_lower()-x.nnz_diag()==0,
1548  "Sparsity error in triu2symm. Found below-diagonal entries in argument: " + x.dim());
1549  return x + x.T() - diag(diag(x));
1550  }

◆ zeros() [1/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::zeros ( casadi_int  nrow = 1,
casadi_int  ncol = 1 
)
inlinestatic

◆ zeros() [2/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::zeros ( const Sparsity sp)
inlinestatic

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

Definition at line 1265 of file generic_matrix.hpp.

1265 { return MatType(sp, 0, false);}

◆ zeros() [3/3]

template<typename MatType >
static MatType casadi::GenericMatrix< MatType >::zeros ( const std::pair< casadi_int, casadi_int > &  rc)
inlinestatic

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

Definition at line 1266 of file generic_matrix.hpp.

1266  {
1267  return zeros(rc.first, rc.second);
1268  }

References casadi::GenericMatrix< MatType >::zeros().


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