---
title: "Matrix, vector, cube and field classes"
output: rmarkdown::html_vignette
bibliography: "references.bib"
vignette: >
  %\VignetteIndexEntry{Matrix, vector, cube and field classes}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

This vignette is adapted from the official Armadillo
[documentation](https://arma.sourceforge.net/docs.html).

# Dense matrix class

`Mat<type>`, `mat` and `cx_mat` are classes for dense matrices, with elements
stored in [column-major ordering](https://en.wikipedia.org/wiki/Column_major)
(e.g., column by column).
  
The root matrix class is `Mat<type>`, where `type` is one of:

  * `float`
  * `double`
  * `std::complex<float>`
  * `std::complex<double>`
  * `short`
  * `int`
  * `long`
  * `unsigned short`
  * `unsigned int`
  * `unsigned long`

For convenience the following typedefs have been defined:
    
* `mat = Mat<double>`
* `dmat = Mat<double>`
* `fmat = Mat<float>`
* `cx_mat = Mat<cx_double>` (`cx_double` is a shortcut for
  `std::complex<double>`)
* `cx_dmat = Mat<cx_double>`
* `cx_fmat = Mat<cx_float>` (`cx_float` is a shortcut for `std::complex<float>`)
* `umat = Mat<uword>` (`uword` is a shortcut for `unsigned int`)
* `imat = Mat<sword>` (`sword` is a shortcut for `signed int`)
    
The `mat` type is used for convenience, and it is possible to use other matrix
types (e.g, `fmat`, `cx_mat`) instead.
  
Matrix types with integer elements (such as `umat` and `imat`) cannot hold
special values such as NaN and Inf.
  
Functions which use LAPACK (generally matrix decompositions) are only valid for
the following matrix types: `mat`, `dmat`, `fmat`, `cx_mat`, `cx_dmat`,
`cx_fmat`.
  
## Constructors
  
  * `mat()`
  * `mat(n_rows, n_cols)`
  * `mat(n_rows, n_cols, fill_form)` (elements are initialised according to
    `fill_form`)
  * `mat(size(X))`
  * `mat(size(X), fill_form)` (elements are initialised according to
    `fill_form`)
  * `mat(mat)`
  * `mat(vec)`
  * `mat(rowvec)`
  * `mat(initializer_list)`
  * `mat(string)`
  * `mat(std::vector)` (treated as a column vector)
  * `mat(sp_mat)` (for converting a sparse matrix to a dense matrix)
  * `cx_mat(mat,mat)` (for constructing a complex matrix out of two real
    matrices)
    
The elements can be explicitly initialised during construction by specifying
`fill_form`, which is one of:
    
  * `fill::zeros` set all elements to 0 (default in cpp11armadillo)
  * `fill::ones` set all elements to 1
  * `fill::eye` set the elements on the main diagonal to 1 and off-diagonal
    elements to 0
  * `fill::randu` set all elements to random values from a uniform distribution
    in the [0,1] interval
  * `fill::randn` set all elements to random values from a normal distribution
    with zero mean and unit variance
  * `fill::value(scalar)` set all elements to specified scalar
  * `fill::none` do not initialise the elements (matrix may have garbage values)
      
For the `mat(string)` constructor, the format is elements separated by spaces,
and rows denoted by semicolons. For example, the 2x2 identity matrix can be
created using `"1 0; 0 1"`. Note that string based initialisation is slower than
directly setting the elements or using element initialisation.
  
Each instance of `mat` automatically allocates and releases internal memory. All
internally allocated memory used by an instance of mat is automatically released
as soon as the instance goes out of scope. For example, if an instance of mat is
declared inside a function, it will be automatically destroyed at the end of the
function. To forcefully release memory at any point, use `.reset()`. Note that
in normal use this is not required.
  
## Advanced constructors
      
* `mat(ptr_aux_mem, n_rows, n_cols, copy_aux_mem = true, strict = false)`
  * Create a matrix using data from writable auxiliary (external) memory, where
    `ptr_aux_mem` is a pointer to the memory. By default the matrix allocates
    its own memory and copies data from the auxiliary memory (for safety).
    However, if `copy_aux_mem` is set to `false`, the matrix will instead
    directly use the auxiliary memory (e.g., no copying). This is faster, but
    can be dangerous unless you know what you are doing.    
  * The `strict` parameter comes into effect only when copy_aux_mem is set to
    `false` (e.g., the matrix is directly using auxiliary memory).
    * When `strict` is set to `false`, the matrix will use the auxiliary memory
      until a size change or an aliasing event.
    * When `strict` is set to `true`, the matrix will be bound to the auxiliary
      memory for its lifetime. The number of elements in the matrix cannot be
      changed.
* `mat(const ptr_aux_mem, n_rows, n_cols)`  
  * Create a matrix by copying data from read-only auxiliary memory, where
    `ptr_aux_mem` is a pointer to the memory
* `mat::fixed<n_rows, n_cols>`
  * Create a fixed size matrix, with the size specified via template arguments.
    Memory for the matrix is reserved at compile time. This is generally faster
    than dynamic memory allocation, but the size of the matrix cannot be changed
    afterwards (directly or indirectly).  
  * For convenience, there are several pre-defined typedefs for each matrix type
    (where the types are: `umat`, `imat`, `fmat`, `mat`, `cx_fmat`, `cx_mat`).
    The typedefs specify a square matrix size, ranging from 2x2 to 9x9. The
    typedefs were defined by appending a two digit form of the size to the
    matrix type. Examples: `mat33` is equivalent to `mat::fixed<3,3>`, and
    `cx_mat44` is equivalent to `cx_mat::fixed<4,4>`.
* `mat::fixed<n_rows, n_cols>(fill_form)`  
  * Create a fixed size matrix, with the elements explicitly initialised
    according to `fill_form`.
* `mat::fixed<n_rows, n_cols>(const ptr_aux_mem)`
  * Create a fixed size matrix, with the size specified via template arguments;
    data is copied from auxiliary memory, where ptr_aux_mem is a pointer to the
    memory.
    
## Examples

```r
set.seed(123)
a <- matrix(runif(25), nrow = 5, ncol = 5)
```

```cpp
[[cpp11::register]] doubles_matrix<> matrix1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a); // convert from R to C++

  double x = A(0, 0);  // access an element on row 1, column 1
  A = A + x; // scalar addition
      
  mat B = A + A; // matrix addition
  mat C = A * B; // matrix multiplication
  mat D = A % B; // element-wise matrix multiplication

  mat res = B + C + D;

  return as_doubles_matrix(res);  // convert from C++ to R
}
      
[[cpp11::register]] list matrix2_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);
  mat B = A + A;

  cx_mat X(A,B); // construct a complex matrix out of two real matrices

  B.zeros();                       // set all elements to zero
  B.set_size(A.n_rows, A.n_cols);  // resize the matrix
  B.ones(5, 6);                    // same as mat B(5, 6, fill::ones)
        
  mat::fixed<5,6> F; // fixed size matrix
      
  double aux_mem[24]; // auxiliary memory
  mat H(&aux_mem[0], 4, 6, false);  // use auxiliary memory

  X = X + F.submat(0, 0, 4, 4) + H(1, 2)

  Mat<double> res_real = real(X);
  Mat<double> res_imag = imag(X);

  writable::list res;
  res.push_back({"real"_nm = as_doubles_matrix(res_real)});
  res.push_back({"imag"_nm = as_doubles_matrix(res_imag)});

  return res;
}
```

# Column vector class

`Col<type>`, `vec` and `cx_vec` are classes for column vectors (dense matrices
with one column).
  
The `Col<type>` class is derived from the `Mat<type>` class and inherits most of
the member functions.
  
For convenience the following typedefs have been defined:
    
* `vec = colvec = Col<double>`
* `dvec = dcolvec = Col<double>`
* `fvec = fcolvec = Col<float>`
* `cx_vec = cx_colvec = Col<[cx_double](#cx_double)>`
* `cx_dvec = cx_dcolvec = Col<[cx_double](#cx_double)>`
* `cx_fvec = cx_fcolvec = Col<[cx_float](#cx_double)>`
* `uvec = ucolvec = Col<[uword](#uword)>`
* `ivec = icolvec = Col<[sword](#uword)>`

The `vec` and `colvec` types have the same meaning and are used interchangeably.
  
The types `vec` or `colvec` are used for convenience. It is possible to use
other column vector types instead (e.g., `fvec` or `fcolvec`).
  
Functions which take `mat` as input can generally also take `Col` as input.
Main exceptions are functions that require square matrices.
  
## Constructors
    
* `vec()`
* `vec(_n_elem_)`
* `vec(_n_elem, fill_form)` (elements are initialised according to `fill_form`)
* `vec(size(X))`
* `vec(size(X), fill_form)` (elements are initialised according to `fill_form`)
* `vec(vec)`
* `vec(mat)` (`std::logic_error` exception is thrown if the given matrix has more than one column)
* `vec(initializer_list)`
* `vec(string)` (elements separated by spaces)
* `vec(std::vector)`
* `cx_vec(vec,vec)` (for constructing a complex vector out of two real vectors)
      
## Advanced constructors
      
* `vec(ptr_aux_mem, number_of_elements, copy_aux_mem = true, strict = false)`:
  * Create a column vector using data from writable auxiliary (external) memory,
    where ptr_aux_mem is a pointer to the memory. By default the vector
    allocates its own memory and copies data from the auxiliary memory (for
    safety). However, if copy_aux_mem is set to `false`, the vector will instead
    directly use the auxiliary memory (e.g., no copying). This is faster, but can
    be dangerous unless you know what you are doing.
  * The `strict` parameter comes into effect only when `copy_aux_mem` is set to
    `false` (e.g., the vector is directly using auxiliary memory).
  * When `strict` is set to `false`, the vector will use the auxiliary memory
    until a size change or an aliasing event.
  * When `strict` is set to `true`, the vector will be bound to the auxiliary
    memory for its lifetime. The number of elements in the vector cannot be
    changed.
* `vec(const ptr_aux_mem, number_of_elements)`: Create a column vector by
    copying data from read-only auxiliary memory, where ptr_aux_mem is a pointer
    to the memory.
* `vec::fixed<number_of_elements>`:
  * Create a fixed size column vector, with the size specified via the template
    argument. Memory for the vector is reserved at compile time. This is
    generally faster than dynamic memory allocation, but the size of the vector
    cannot be changed afterwards (directly or indirectly).  
  * For convenience, there are several pre-defined typedefs for each vector type
    (where the types are: `uvec`, `ivec`, `fvec`, `vec`, `cx_fvec`, `cx_vec` as
    well as the corresponding `colvec` versions). The pre-defined typedefs
    specify vector sizes ranging from 2 to 9. The typedefs were defined by
    appending a single digit form of the size to the vector type.
    Examples: `vec3` is equivalent to `vec::fixed<3>`, and `cx_vec4` is
    equivalent to `cx_vec::fixed<4>`.
* `vec::fixed<number_of_elements>(fill_form)`: Create a fixed size column
    vector, with the elements explicitly initialised according to `fill_form`.
* `vec::fixed<number_of_elements>(const ptr_aux_mem)`: Create a fixed size
    column vector, with the size specified via the template argument. The data
    is copied from auxiliary memory, where `ptr_aux_mem` is a pointer to the
    memory.
      
## Examples
  
```r
set.seed(123)
x <- runif(10)
y <- rep(1, 10)
```

```cpp
[[cpp11::register]] doubles column1_(const doubles& x, const doubles& y) {
  vec X = as_Col(x); // convert from R to C++
  vec Y = as_Col(y);

  mat A(10, 10, fill::randu);
  vec Z = A.col(5); // extract a column vector

  Z = Z + Y + X;

  return as_doubles(Z); // convert from C++ to R
}
```

# Row vector class

`Row<type>`, `rowvec` and `cx_rowvec` are classes for row vectors (dense
matrices with one row).
  
The template `Row<type>` class is derived from the `Mat<type>` class and
inherits most of the member functions.
  
For convenience the following typedefs have been defined:

* `rowvec = Row<double>`
* `drowvec = Row<double>`
* `frowvec = Row<float>`
* `cx_rowvec = Row<cx_double>`
* `cx_drowvec = Row<cx_double>`
* `cx_frowvec = Row<cx_float>`
* `urowvec = Row<uword>`
* `irowvec = Row<sword>`
    
The `rowvec` type is used for convenience. It is possible to use other row
vector types instead (e.g., `frowvec`).
  
Functions which take `mat` as input can generally also take `Row` as input. Main
exceptions are functions which require square matrices.
  
## Constructors
    
* `rowvec()`
* `rowvec(n_elem)`
* `rowvec(n_elem, fill_form)` (elements are initialised according to `fill_form`)
* `rowvec(size(X))`
* `rowvec(size(X), fill_form)` (elements are initialised according to `fill_form`)
* `rowvec(rowvec)`
* `rowvec(mat)` (`std::logic_error` exception is thrown if the given matrix has more than one row)
* `rowvec(initializer_list)`
* `rowvec(string)` (elements separated by spaces)
* `rowvec(std::vector)`
* `cx_rowvec(rowvec,rowvec)` (for constructing a complex row vector out of two real row vectors)
    
## Advanced constructors
      
* `rowvec(ptr_aux_mem, number_of_elements, copy_aux_mem = true, strict = false)`  
  * Create a row vector using data from writable auxiliary (external) memory,
    where ptr_aux_mem is a pointer to the memory. By default the vector
    allocates its own memory and copies data from the auxiliary memory (for
    safety). However, if copy_aux_mem is set to `false`, the vector will instead
    directly use the auxiliary memory (e.g., no copying); this is faster, but can
    be dangerous unless you know what you are doing. 
  * The `strict` parameter comes into effect only when copy_aux_mem is set to
    `false` (e.g., the vector is directly using auxiliary memory).
    * When `strict` is set to `false`, the vector will use the auxiliary memory
      until a size change or an aliasing event.
    * When `strict` is set to `true`, the vector will be bound to the auxiliary
      memory for its lifetime. The number of elements in the vector cannot be
      changed.
* `rowvec(const ptr_aux_mem, number_of_elements)`  
  * Create a row vector by copying data from read-only auxiliary memory, where
    `ptr_aux_mem` is a pointer to the memory.
* `rowvec::fixed<number_of_elements>`
  * Create a fixed size row vector, with the size specified via the template
    argument. Memory for the vector is reserved at compile time. This is
    generally faster than dynamic memory allocation, but the size of the vector
    cannot be changed afterwards (directly or indirectly).
  * For convenience, there are several pre-defined typedefs for each vector type
    (where the types are: `urowvec`, `irowvec`, `frowvec`, `rowvec`,
    `cx_frowvec`, `cx_rowvec`). The pre-defined typedefs specify vector sizes
    ranging from 2 to 9. The typedefs were defined by appending a single digit
    form of the size to the vector type. Examples: `rowvec3` is equivalent to
    `rowvec::fixed<3>`, and `cx_rowvec4` is equivalent to `cx_rowvec::fixed<4>`.
* `rowvec::fixed<number_of_elements>(fill_form)`
  * Create a fixed size row vector, with the elements explicitly initialised
    according to `fill_form`.
* `rowvec::fixed<number_of_elements>(const ptr_aux_mem)`  
  * Create a fixed size row vector, with the size specified via the template
    argument. The data is copied from auxiliary memory, where `ptr_aux_mem` is a
    pointer to the memory.
    
## Examples
    
⚠️Important⚠️: 'cpp11armadillo' is an opinionated package and it follows the
notation from Econometrics by Bruce E. Hansen. It intentionally exports/imports
matrices and column vectors. You can use row vectors in the functions, but
the communication between R and C++ does not accept row vectors unless you
transpose or convert those to matrices.

```r
set.seed(123)
x <- runif(10)
y <- rep(1, 10)
```

```cpp
[[cpp11::register]] doubles row1_(const doubles& x, const doubles& y) {
  vec X = as_Col(x);  // convert from R to C++
  vec Y = as_Col(y);

  mat A(10, 10, fill::randu);
  
  rowvec Z = A.row(5);  // extract a row vector
  Z = Z + Y.t() + X.t(); // transpose Y and X to be able to sum

  vec res = Z.t();
  
  return as_doubles(res);  // convert from C++ to R
}
```

# Cube class

`Cube<type>`, `cube` and `cx_cube` are classes for cubes, also known as quasi
3rd order tensors or "3D matrices".
  
The data is stored as a set of slices (matrices) stored contiguously within
memory. Within each slice, elements are stored with column-major ordering (e.g.,
column by column)
  
The root cube class is `Cube<type>`, where `type` is one of: `float`, `double`,
`std::complex<float>`, `std::complex<double>`, `short`, `int`, `long` and
unsigned versions of `short`, `int`, `long`.
  
For convenience the following typedefs have been defined:
    
* `cube = Cube<double>`
* `dcube = Cube<double>`
* `fcube = Cube<float>`
* `cx_cube = Cube<cx_double>`
* `cx_dcube = Cube<cx_double>`
* `cx_fcube = Cube<cx_float>`
* `ucube = Cube<uword>`
* `icube = Cube<sword>`
    
The `cube` type is used for convenience. It is possible to use other types
instead (e.g., `fcube`).

Each cube slice can be interpreted as a matrix, hence functions which take
`Mat` as input can generally also take cube slices as input.
  
## Constructors

* `cube()`
* `cube(n_rows, n_cols, n_slices_)`
* `cube(n_rows, n_cols, n_slices, fill_form)` (elements are initialised according to `fill_form`)
* `cube(size(X))`
* `cube(size(X), fill_form)` (elements are initialised according to `fill_form`)
* `cube(cube)`
* `cx_cube(cube, cube)` (for constructing a complex cube out of two real cubes)
      
The elements can be explicitly initialised during construction by specifying `fill_form`, which is one of:
    
* `fill::zeros`: set all elements to 0 (default in cpp11armadillo)
* `fill::ones`: set all elements to 1
* `fill::randu`: set all elements to random values from a uniform distribution in the [0,1] interval
* `fill::randn`: set all elements to random values from a normal  distribution with zero mean and unit variance
* `fill::value(scalar)`: set all elements to specified scalar
* `fill::none`: do not initialise the elements (cube may have garbage values)
    
Each instance of `cube` automatically allocates and releases internal memory.
All internally allocated memory used by an instance of `cube` is automatically
released as soon as the instance goes out of scope. For example, if an instance
of `cube` is declared inside a function, it will be automatically destroyed at
the end of the function. To forcefully release memory at any point, use
`.reset()` note that in normal use this is not required.
  
## Advanced constructors
      
* `cube::fixed<n_rows, n_cols, n_slices>`: Create a fixed size cube, with the
    size specified via template arguments. Memory for the cube is reserved at
    compile time. This is generally faster than dynamic memory allocation, but
    the size of the cube cannot be changed afterwards (directly or indirectly).
* `cube(ptr_aux_mem, n_rows, n_cols, n_slices, copy_aux_mem = true, strict = false)`: 
  * Create a cube using data from writable auxiliary (external) memory, where
    `ptr_aux_mem` is a pointer to the memory. By default the cube allocates its
    own memory and copies data from the auxiliary memory (for safety). However,
    if `copy_aux_mem` is set to `false`, the cube will instead directly use the
    auxiliary memory (e.g., no copying). This is faster, but can be dangerous
    unless you know what you are doing.
  * The `strict` parameter comes into effect only when `copy_aux_mem` is set to
    `false` (e.g., the cube is directly using auxiliary memory).
  * When `strict` is set to `false`, the cube will use the auxiliary memory
    until a size change or an aliasing event.
  * When `strict` is set to `true`, the cube will be bound to the auxiliary
    memory for its lifetime. The number of elements in the cube cannot be
    changed.
* `cube(const ptr_aux_mem, n_rows, n_cols, n_slices)`: Create a cube by copying
    data from read-only auxiliary memory, where `ptr_aux_mem` is a pointer to
    the memory.
    
## Examples

```r
set.seed(123)
a <- matrix(runif(4), nrow = 2, ncol = 2)
b <- matrix(rnorm(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> cube1_(const doubles_matrix<>& a,
                                                const doubles_matrix<>& b) {
  mat A = as_Mat(a);  // convert from R to C++
  mat B = as_Mat(b);

  cube X(A.n_rows, A.n_cols, 2);  // create a cube with 2 slices
  X.slice(0) = A;                 // copy A into first slice
  X.slice(1) = B;                 // copy B into second slice

  cube Y = X + X;  // cube addition
  cube Z = X % X;  // element-wise cube multiplication

  mat res = Y.slice(0) + Z.slice(1);

  return as_doubles_matrix(res);  // convert from C++ to R
}
```
  
The size of individual slices cannot be changed. The following will not work:
        
```cpp
cube c(5,6,7);
c.slice(0) = randu<mat>(10,20); // wrong size
```

# Field class

`field<object_type>` is a class for storing arbitrary objects in matrix-like or
cube-like layouts.
  
It is similar to a matrix or cube, but instead of each element being a scalar,
each element can be a vector, or matrix, or cube. This is similar to a list in
R.
  
Each element can have an arbitrary size (e.g., in a field of matrices, each
matrix can have a unique size).
  
## Constructors

`object_type` is another class (e.g., `vec`, `mat`, `std::string`, etc)
    
* `field<object_type>()`
* `field<object_type>(n_elem)`
* `field<object_type>(n_rows, n_cols)`  
* `field<object_type>(n_rows, n_cols, n_slices)`
* `field<object_type>(size(X))`
* `field<object_type>(field<object_type>)`

## Examples

```r
set.seed(123)
a <- matrix(runif(4), nrow = 2, ncol = 2)
b <- matrix(rnorm(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> field1_(const doubles_matrix<>& a,
                                                 const doubles_matrix<>& b) {
  mat A = as_Mat(a);  // convert from R to C++
  mat B = as_Mat(b);

  field<mat> F(A.n_rows, A.n_cols, 3);  // create a field with 2 matrices
  F(0) = A;                             // copy A into first location
  F(1) = B;                             // copy B into second location
  F(2) = F(0) + F(1);                   // matrix addition

  mat res = F(0) + F(1) + F(2).t();

  return as_doubles_matrix(res);  // convert from C++ to R
}
```

## Caveat

To store a set of matrices of the same size, the `Cube` class is more efficient.

# Member Functions and Variables

| Function/Variable | Description                                                                 |
|-------------------|-----------------------------------------------------------------------------|
| `.n_rows`         | number of rows                                                              |
| `.n_cols`         | number of columns                                                           |
| `.n_elem`         | number of elements                                                          |
| `.n_slices`       | number of slices                                                            |
| `()`              | element access                                                              |
| `[]`              | element access                                                              |
| `.at()`           | element access                                                              |
| `.zeros`          | set all elements to zero                                                    |
| `.ones`           | set all elements to one                                                     |
| `.eye`            | set elements along main diagonal to one and off-diagonal elements to zero   |
| `.randu`          | set all elements to random values from a uniform distribution               |
| `.randn`          | set all elements to random values from a normal distribution                |
| `.fill`           | set all elements to specified value                                         |
| `.imbue`          | imbue (fill) with values provided by functor or lambda function             |
| `.clean`          | replace elements below a threshold with zeros                               |
| `.replace`        | replace specific elements with a new value                                  |
| `.clamp`          | clamp values to lower and upper limits                                      |
| `.transform`      | transform each element via functor or lambda function                       |
| `.for_each`       | apply a functor or lambda function to each element                          |
| `.set_size`       | change size without keeping elements (fast)                                 |
| `.reshape`        | change size while keeping elements                                          |
| `.resize`         | change size while keeping elements and preserving layout                    |
| `.copy_size`      | change size to be same as given object                                      |
| `.reset`          | change size to empty                                                        |
| `.diag`           | read/write access to matrix diagonals                                       |
| `.each_col`       | vector operations applied to each column of matrix (aka "broadcasting")     |
| `.each_row`       | vector operations applied to each row of matrix (aka "broadcasting")        |
| `.each_slice`     | matrix operations applied to each slice of cube (aka "broadcasting")        |
| `.set_imag`       | set imaginary part                                                          |
| `.set_real`       | set real part                                                               |
| `.insert_rows`    | insert vector/matrix/cube at specified row                                  |
| `.insert_cols`    | insert vector/matrix/cube at specified column                               |
| `.insert_slices`  | insert vector/matrix/cube at specified slice                                |
| `.shed_rows`      | remove specified rows                                                       |
| `.shed_cols`      | remove specified columns                                                    |
| `.shed_slices`    | remove specified slices                                                     |
| `.swap_rows`      | swap specified rows                                                         |
| `.swap_cols`      | swap specified columns                                                      |
| `.swap`           | swap contents with given object                                             |
| `.memptr`         | raw pointer to memory                                                       |
| `.colptr`         | raw pointer to memory used by specified column                              |
| `.as_col`         | return flattened matrix column as column vector                             |
| `.as_row`         | return flattened matrix row as row vector                                   |
| `.col_as_mat`     | return matrix representation of cube column                                 |
| `.row_as_mat`     | return matrix representation of cube row                                    |
| `.as_dense`       | return dense vector/matrix representation of sparse matrix expression       |
| `.t`              | return matrix transpose                                                     |
| `.st`             | return matrix conjugate transpose                                           |
| `.i`              | return inverse of square matrix                                             |
| `.min`            | return minimum value                                                        |
| `.max`            | return maximum value                                                        |
| `.index_min`      | return index of minimum value                                               |
| `.index_max`      | return index of maximum value                                               |
| `.eval`           | force evaluation of delayed expression                                      |
| `.in_range`       | check whether given location or span is valid                               |
| `.is_empty`       | check whether object is empty                                               |
| `.is_vec`         | check whether matrix is a vector                                            |
| `.is_sorted`      | check whether vector or matrix is sorted                                    |
| `.is_trimatu`     | check whether matrix is upper triangular                                    |
| `.is_trimatl`     | check whether matrix is lower triangular                                    |
| `.is_diagmat`     | check whether matrix is diagonal                                            |
| `.is_square`      | check whether matrix is square sized                                        |
| `.is_symmetric`   | check whether matrix is symmetric                                           |
| `.is_hermitian`   | check whether matrix is hermitian                                           |
| `.is_sympd`       | check whether matrix is symmetric/hermitian positive definite               |
| `.is_zero`        | check whether all elements are zero                                         |
| `.is_finite`      | check whether all elements are finite                                       |
| `.has_inf`        | check whether any element is +/- infinity                                   |
| `.has_nan`        | check whether any element is NaN                                            |

# Attributes

`n_*` provides information for different objects:

* `.n_rows` number of rows for `Mat`, `Col`, `Row`, `Cube`, `field`, and
  `SpMat`.
* `.n_cols` number of columns for `Mat`, `Col`, `Row`, `Cube`, `field`, and
  `SpMat`.
* `.n_elem` total number of elements for `Mat`, `Col`, `Row`, `Cube`, `field`,
  and `SpMat`.
* `.n_slices` number of slices for `Cube` and `field`.
* `.n_nonzero` number of non-zero elements for `SpMat`.
  
For the `Col` and `Row` classes, `n_elem` also indicates vector length.

The variables are read-only and of type `uword`. To change the size, use
`set_size`, `copy_size`, `zeros_member`, `ones_member`, or `reset`.

To avoid compiler warnings about implicit conversion when operating `uword` with
`integers`/`doubles` to pass data to R, converte `uword` to `int` with
`static_cast<int>` or declare these as `int`.

## Examples

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] integers attr1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++

  // uword or int can be used
  int n_rows = A.n_rows;  // number of rows
  int n_cols = A.n_cols;  // number of columns
  int n_elem = A.n_elem;  // number of elements

  writable::integers res({n_rows, n_cols, n_elem});
  res.attr("names") = strings({"n_rows", "n_cols", "n_elem"});

  return res;
}
```

# Element/object access

Provide access to individual elements or objects stored in a container object
(e.g., `Mat`, `Col`, `Row`, `Cube`, `field`).
      
* `(i)` For `vec` and `rowvec`, access the element stored at index `i`. For
  `Mat`, `Cube` and `field`, access the element/object stored at index `i` under
  the assumption of a flat layout, with column-major ordering of data (e.g.,
  column by column). An exception is thrown if the requested element is out of
  bounds.
* `.at(i)` or `[i]` As for `(i)`, but without a bounds check. Not recommended.
* `(r,c)` For `Mat` and 2D field classes, access the element/object stored at 
  row `r` and column `c`. An exception is thrown if the requested element is out
  of bounds.
* `.at(r,c)` As for `(r,c)`, but without a bounds check. Not recommended.
* `(r,c,s)` For `Cube` and 3D field classes, access the element/object stored at
  row `r`, column `c`, and slice `s`. An exception is thrown if the requested
  element is out of bounds.
* `.at(r,c,s)` As for `(r,c,s)`, but without a bounds check. Not recommended.

## Examples

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> access1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++
  A(1,1) = 123.0;     // set element at row 2, column 2

  vec B(2, fill::randu);

  double x = A(0,1);  // copy element at row 1, column 2 to a double
  double y = B(1);    // copy element at coordinate 2 to a double

  uword i, j; // int also works
  uword N = A.n_rows;
  uword M = A.n_cols;

  for(i = 0; i < N; ++i) {
    for(j = 0; j < M; ++j) {
      A(i,j) = A(i,j) + x + y;
    }
  }

  return as_doubles_matrix(A);  // convert from C++ to R
}
```

## Caveats

For `.at()` or `[i]`, `.at(r,c)` and `.at(r,c,s)`:

  * Indexing in C++ starts at 0
  * Accessing elements without bounds checks is slightly faster, but is not
    recommended until your code has been thoroughly debugged first
  * Accessing elements via `[r,c]` and `[r,c,s]` does not work correctly in C++;
    instead use `(r,c)` and `(r,c,s)`

The indices of elements are specified via the `uword` type, which is a `typedef`
for an unsigned integer type. When using loops to access elements, it more
efficient to use `uword` instead of `int`.

# Element initialisation

Set elements in `Mat`, `Col` and `Row` via braced initialiser lists.
  
## Examples
    
```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> initialization1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++
  mat B = {{1, 2}, {3, 4}}; // create new matrix
  vec C = {1, 2}; // create new column vector

  // sum C to the diagonal of A
  A(0,0) = A(0,0) + C(0);
  A(1,1) = A(1,1) + C(1);

  mat D = A + B;
  
  return as_doubles_matrix(D);  // convert from C++ to R
}
```

# Zeros

Set the elements of an object to zero, optionally first changing the size to
specified dimensions.

`.zeros()` (member function of `Mat`, `Col`, `Row`, `SpMat`, `Cube`)
`.zeros(n_elem)` (member function of `Col` and `Row`)
`.zeros(n_rows, n_cols)` (member function of `Mat` and `SpMat`)
`.zeros(n_rows, n_cols, n_slices)` (member function of `Cube`)
`.zeros(size(X))` (member function of `Mat`, `Col`, `Row`, `Cube`, `SpMat`)
  
## Examples
    
```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> zeros1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++
  A.zeros();          // set all elements to zero

  mat B;
  B.zeros(size(A)); // set size to be the same as A and set all elements to zero

  mat C(A.n_rows, A.n_cols, fill::zeros);

  mat D = A + B + C;

  return as_doubles_matrix(D);  // convert from C++ to R
}
```

# Ones

Set all the elements of an object to one, optionally first changing the size to
specified dimensions.

| Function                                 | Mat | Col | Row | Cube |
|------------------------------------------|-----|-----|-----|------|
| `.ones()`                                | ✓   | ✓   | ✓   | ✓    |
| `.ones(n_elem)`                          |     | ✓   | ✓   |      |
| `.ones(n_rows, n_cols)`                  | ✓   |     |     |      |
| `.ones(n_rows, n_cols, n_slices)`        |     |     |     | ✓    |
| `.ones(size(X))`                         | ✓   | ✓   | ✓   | ✓    |
  
## Examples
    
```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> ones1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++
  A.ones();          // set all elements to zero

  mat B;
  B.ones(size(A)); // set size to be the same as A and set all elements to zero

  mat C(A.n_rows, A.n_cols, fill::ones);

  mat D = A + B + C;

  return as_doubles_matrix(D);  // convert from C++ to R
}
```

# Eye

`.eye()` is member function of `Mat` and `SpMat`. `.eye(n_rows, n_cols)` sets
the elements along the main diagonal to one and off-diagonal elements to zero,
optionally first changing the size to specified dimensions. `.eye(size(X))`
creates an identity matrix is generated when `n_rows = n_cols`.
  
## Examples
    
```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> eye1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++
  A.eye();            // create an identity matrix

  mat B;
  B.eye(size(A)); // another identity matrix

  uword N = A.n_rows;
  uword M = A.n_cols;
  mat C(N, M, fill::randu);
  C.eye(N, M); // yet another identity matrix

  mat D = A + B + C;

  return as_doubles_matrix(D);  // convert from C++ to R
}
```

# Random uniform

Set all the elements to random values from a uniform distribution in the [0,1]
interval, optionally first changing the size to specified dimensions.

For complex elements, the real and imaginary parts are treated separately.

| Function/Method                          | Mat | Col | Row | Cube |
|------------------------------------------|-----|-----|-----|------|
| `.randu()`                               | ✓   | ✓   | ✓   | ✓    |
| `.randu(n_elem)`                         |     | ✓   | ✓   |      |
| `.randu(n_rows, n_cols)`                 | ✓   |     |     |      |
| `.randu(n_rows, n_cols, n_slices)`       |     |     |     | ✓    |
| `.randu(size(X))`                        | ✓   | ✓   | ✓   | ✓    |

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> randu1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++

  mat B;
  B.randu(size(A)); // random uniform matrix with the same size as A

  mat C(A.n_rows, A.n_cols, fill::randu);

  mat D = A + B + C;

  return as_doubles_matrix(D);  // convert from C++ to R
}

[[cpp11::register]] doubles_matrix<> randu2_(const int& n) {
  GetRNGstate();  // Ensure R's RNG state is synchronized
  mat y(n, n);
  ::arma_rng::randu<double>::fill(y.memptr(), y.n_elem);
  PutRNGstate();

  return as_doubles_matrix(y);
}
```

# Normal distribution

Set all the elements to random values from a normal distribution with zero mean
and unit variance, optionally first changing the size to specified dimensions.

For complex elements, the real and imaginary parts are treated separately.

| Function/Method                          | Mat | Col | Row | Cube |
|------------------------------------------|-----|-----|-----|------|
| `.randn()`                               | ✓   | ✓   | ✓   | ✓    |
| `.randn(n_elem)`                         |     | ✓   | ✓   |      |
| `.randn(n_rows, n_cols)`                 | ✓   |     |     |      |
| `.randn(n_rows, n_cols, n_slices)`       |     |     |     | ✓    |
| `.randn(size(X))`                        | ✓   | ✓   | ✓   | ✓    |

## Examples

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> randn1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++

  mat B;
  B.randn(size(A)); // random normal matrix with the same size as A

  mat C(A.n_rows, A.n_cols, fill::randn);

  mat D = A + B + C;

  return as_doubles_matrix(D);  // convert from C++ to R
}

[[cpp11::register]] doubles_matrix<> randn2_(const int& n) {
  GetRNGstate();  // Ensure R's RNG state is synchronized
  mat y(n, n);
  ::arma_rng::randn<double>::fill(y.memptr(), y.n_elem);
  PutRNGstate();
  
  return as_doubles_matrix(y);
}
```

# Fill

Sets the elements to a specified value

`.fill(value)` is a member function of `Mat`, `Col`, `Row`, `Cube`, `field`.
  
The type of value must match the type of elements used by the container object
(e.g., for `Mat` the type is `double`)
  
## Examples

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> fill1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++

  uword N = A.n_rows;
  uword M = A.n_cols;

  mat B(size(A), fill::value(200.0));  // create a matrix filled with 200.0
  mat C(N, M, fill::value(100.0));     // matrix filled with 100.0
  mat D(N, M, fill::zeros);            // matrix filled with zeros
  mat E(N, M, fill::ones);             // matrix filled with ones

  mat F = A + B + C + D + E;

  return as_doubles_matrix(F);  // convert from C++ to R
}
```

# Imbue

`.imbue(functor)` is a member function of `Mat`, `Col`, `Row` and `Cube`, it
fills the elements with values provided by a functor. The argument can be a
functor or lambda function.

For matrices, filling is done column-by-column (e.g., column 0 is filled, then
column 1, etc.)
  
For cubes, filling is done slice-by-slice, with each slice treated as a matrix
  
## Examples

```r
a <- matrix(runif(4), nrow = 2, ncol = 2)
```

```cpp
[[cpp11::register]] doubles_matrix<> imbue1_(const doubles_matrix<>& a) {
  mat A = as_Mat(a);  // convert from R to C++

  std::mt19937 engine;  // Mersenne twister random number engine
  std::uniform_real_distribution<double> distr(0.0, 1.0);

  mat B(size(A), fill::none);                // create an empty matrix
  B.imbue([&]() { return distr(engine); });  // fill with random values

  mat C = A + B;

  return as_doubles_matrix(C);  // convert from C++ to R
}

[[cpp11::register]] doubles_matrix<> imbue2_(const doubles_matrix<>& a) {
  GetRNGstate();  // Ensure R's RNG state is synchronized

  mat A = as_Mat(a);  // Convert from R to C++

  mat B(size(A), fill::none);  // Create an empty matrix
  B.imbue([]() { return unif_rand(); });  // Fill with random values

  mat C = A + B;

  PutRNGstate();

  return as_doubles_matrix(C);  // Convert from C++ to R
}
```

# Clean

`.clean(threshold)` is a member function of `Mat`, `Col`, `Row`, `Cube`, and
`SpMat`. It can be used to sparsify a matrix, in the sense of zeroing values
with small magnitudes.

* For objects with non-complex elements: each element with an absolute value less
or equal to the threshold is replaced by zero.
* For objects with complex elements: for each element, each component (real and
imaginary) with an absolute value less or equal to the threshold is replaced by
zero.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> clean1_(const int& n) {
  mat A(n, n, fill::randu); // create a random matrix

  A(0, 0) = datum::eps; // set the diagonal with small values (+/- epsilon)
  A(1, 1) = -datum::eps;

  A.clean(datum::eps); // set elements with small values to zero

  return as_doubles_matrix(A); // Convert from C++ to R
}
```

## Caveat

To explicitly convert from dense storage to sparse storage, use the `SpMat`.

# Replace

`.replace( old_value, new_value )` is a member function of `Mat`, `Col`, `Row`,
`Cube`, and `SpMat`.
  
For all elements equal to `old_value`, set them to `new_value`.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> replace1_(const int& n) {
  mat A(n, n, fill::randu); // create a random matrix

  A.diag().fill(datum::nan); // set the diagonal with NaN values
  A.replace(datum::nan, 0);  // replace each NaN with 0

  return as_doubles_matrix(A); // Convert from C++ to R
}
```

## Caveats

* The type of `old_value` and `new_value` must match the type of elements used
  by the container object (e.g., for `Mat` the type is `double`).
* Floating point numbers (`float` and `double`) are approximations due to their
  [limited precision](https://en.wikipedia.org/wiki/Floating-point_arithmetic).
* For sparse matrices (`SpMat`), replacement is not done when `old_value = 0`.

# Clamp

`.clamp(min_value, max_value)` is a member function of `Mat`, `Col`, `Row`,
`Cube` and `SpMat` that transforms all values lower than `min_val` to `min_val`,
and all values higher than `max_val` to `max_val`.
  
* For complex elements, the real and imaginary components are clamped
  separately.
* For sparse matrices, clamping is applied only to the non-zero elements.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> clamp1_(const int& n) {
  mat A(n, n, fill::randu); // create a random matrix
  A.diag().fill(0.1);       // set the diagonal with 0.1 values

  A.clamp(0.2, 0.8); // clamp values to the [0.2, 0.8] interval

  return as_doubles_matrix(A); // Convert from C++ to R
}
```

# Transform

`.transform(functor)` is a member function of `Mat`, `Col`, `Row`, `Cube`, and
`SpMat`. The argument can be a functor or lambda function.

* For dense matrices, transformation is done column-by-column for all elements.
* For sparse matrices, transformation is done column-by-column for non-zero
  elements.
* For cubes, transformation is done slice-by-slice, with each slice treated as a
  matrix.
  
## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> transform1_(const int& n) {
  mat A(n, n, fill::ones);  // create a matrix filled with ones
  A.transform([](double val) { return (val + 122.0); });
  return as_doubles_matrix(A); // Convert from C++ to R
}
```

# For each

`.for_each(functor)` is a member function of `Mat`, `Col`, `Row`, `Cube`,
`SpMat`, and `field`. The argument can be a functor or lambda function.

* For dense matrices and fields, the processing is done column-by-column for all
  elements.
* For sparse matrices, the processing is done column-by-column for non-zero
  elements.
* For cubes, processing is done slice-by-slice, with each slice treated as a
  matrix.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> for_each1_(const int& n) {
  // add 122 to each element in a dense matrix, the '&' is important
  mat D(n, n, fill::ones);
  D.for_each([](mat::elem_type& val) { val += 122.0; });

  // add 122 to each non-zero element in a sparse matrix
  sp_mat S;
  S.sprandu(n, n, 1.0);
  S.for_each([](sp_mat::elem_type& val) { val += 123.0; });

  // set the size of all matrices in a field
  field<mat> F(2, 2);
  F.for_each([n](mat& X) { X.zeros(n, n); });  // capture n for the lambda

  mat res = D + S + F(0) + F(1);

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Set size

Change the size of an object, without explicitly preserving data and without
initialising the elements (e.g., elements may contain garbage values, including
`NaN`).

* `.set_size(n_elem)` (member function of `Col`, `Row`, `field`)
* `.set_size(n_rows, n_cols)` (member function of `Mat`, `SpMat`, `field`)
* `.set_size(n_rows, n_cols, n_slices)` (member function of `Cube` and `field`)
* `.set_size(size(X))` (member function of `Mat`, `Col`, `Row`, `Cube`, `SpMat`,
  `field`)

To initialise the elements to zero while changing the size, use `.zeros()`
instead. To explicitly preserve data while changing the size, use `.reshape()`
or `.resize()` instead.
  
## Examples

```cpp
[[cpp11::register]] doubles set_size1_(const int& n) {
  mat A;
  A.set_size(n, n);  // or:  mat A(n, n, fill::none);

  mat B;
  B.set_size(size(A));  // or:  mat B(size(A), fill::none);

  vec C;
  C.set_size(n);  // or:  vec v(n, fill::none);

  A.fill(1.0);  // set all elements to 1.0
  B.fill(2.0);  // set all elements to 2.0
  C.fill(3.0);  // set all elements to 3.0

  vec res = A.col(0) + B.col(1) + C;

  return as_doubles(res);  // Convert from C++ to R
}
```

# Reshape

Recreate an object according to given size specifications, with the elements
taken from the previous version of the object in a column-wise manner. The
elements in the generated object are placed column-wise (e.g., the first column
is filled up before filling the second column)

* `.reshape(n_rows, n_cols)` (member function of `Mat` and `SpMat`)
* `.reshape(n_rows, n_cols, n_slices)` (member function of `Cube`)
* `.reshape(size(X))` (member function of `Mat`, `Cube`, `SpMat`)

The layout of the elements in the recreated object will be different to the
layout in the previous version of the object
  
If the total number of elements in the previous version of the object is less
than the specified size, the extra elements in the recreated object are set to
zero
  
If the total number of elements in the previous version of the object is greater
than the specified size, only a subset of the elements is taken
  
## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> reshape1_(const int& n) {
  mat A(n + 1, n - 1, fill::randu);
  A.reshape(n - 1, n + 1);
  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

## Caveats

* `.reshape()` is considerably slower than `.set_size()`.
* To change the size without preserving data, use `.set_size()`.
* To grow/shrink the object while preserving the elements and the layout of the
  elements, use `.resize()`
* To flatten a matrix into a vector, use `vectorise()` or
  `.as_col()`/`.as_row()`.

# Resize

Resize an object according to given size specifications, while preserving the
elements and the layout of the elements. It can be used for growing or shrinking
an object (e.g., adding/removing rows, and/or columns, and/or slices).

* `.resize(n_elem)`: member function of `Col`, `Row`.
* `.resize(n_rows, n_cols)`: member function of `Mat` and `SpMat`.
* `.resize(n_rows, n_cols, n_slices)`: member function of `Cube`.
* `.resize(size(X))`: member function of `Mat`, `Col`, `Row`, `Cube`, `SpMat`.

## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> resize1_(const int& n) {
  mat A(n + 1, n - 1, fill::randu);
  A.resize(n - 1, n + 1);
  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

## Caveats

* `.resize()` is considerably slower than `.set_size()`.
* to change the size without preserving data, `.set_size()` instead.

# Copy size

`.copy_size(A)` sets the size of a matrix/vector/cube to be the same as
matrix/vector/cube `A`.

## Examples

```cpp
[[cpp11::register]] integers copy_size1_(const int& n) {
  mat A(n, n, fill::randu);

  mat B;
  B.copy_size(A);

  int N = B.n_rows;
  int M = B.n_cols;
  
  writable::integers res({N, M});
  res.attr("names") = strings({"n_rows", "n_cols"});

  return as_integers(res);  // Convert from C++ to R
}
```

## Caveat

To set the size of an object `B`, `A` must be of the same type as `B`. For
example, the size of a matrix cannot be set by providing a cube.

# Reset

`.reset()` sets a matrix/vector size to zero (the object will have no elements).
  
## Examples

```cpp
[[cpp11::register]] integers reset1_(const int& n) {
  mat A(n, n, fill::randu);
  A.reset();

  int N = A.n_rows;
  int M = A.n_cols;
  
  writable::integers res({N, M});
  res.attr("names") = strings({"n_rows", "n_cols"});

  return as_integers(res);  // Convert from C++ to R
}
```

# Submatrix views

A collection of member functions of `Mat`, `Col` and `Row` classes that provide
read/write access to submatrix views.
      
## Contiguous views for matrix

* `X.col(col_number)`
* `X.row(row_number)`
* `X.cols(first_col, last_col)`  
* `X.rows(first_row, last_row)`
* `X.submat(first_row, first_col, last_row, last_col)`
* `X(span(first_row, last_row), span(first_col, last_col))`
* `X(first_row, first_col, size(n_rows, n_cols))`
* `X(first_row, first_col, size(Y))` (`Y` is a matrix)
* `X(span(first_row, last_row), col_number)`
* `X(row_number, span(first_col, last_col))`
* `X.head_cols(number_of_cols)`
* `X.head_rows(number_of_rows)`
* `X.tail_cols(number_of_cols)`
* `X.tail_rows(number_of_rows)`
* `X.unsafe_col(col_number)` (use with caution)

## Contiguous views for vector

* `Y(span(first_index, last_index))`
* `Y.subvec(first_index, last_index)`
* `Y.subvec(first_index, size(X))` (`X` is a vector)          
* `Y.head(number_of_elements)`
* `Y.tail(number_of_elements)`
    
## Non-contiguous views for matrix or vector:
          
* `X.elem(vector_of_indices)`
* `X(vector_of_indices)`
* `X.cols(vector_of_column_indices)`
* `X.rows(vector_of_row_indices)`
* `X.submat(vector_of_row_indices, vector_of_column_indices)`
* `X(vector_of_row_indices, vector_of_column_indices)`
    
Instances of `span(start, end)` can be replaced by `span::all_` to indicate the
entire range.
  
For functions requiring one or more vector of indices, for example
`X.submat(vector_of_row_indices, vector_of_column_indices)`, each vector of
indices must be of type `uvec`.
  
In the function `X.elem(vector_of_indices)`, elements specified in
`vector_of_indices` are accessed. `X` is interpreted as one long vector, with
column-by-column ordering of the elements of `X`. The `vector_of_indices` must
evaluate to a vector of type `uvec` (e.g., generated by the `find()` function).
The aggregate set of the specified elements is treated as a column vector (e.g.,
the output of `X.elem()` is always a column vector).
  
The function `.unsafe_col()` is provided for speed reasons and should be used
only if you know what you are doing. It creates a seemingly independent `Col`
vector object (e.g., `vec`), but uses memory from the existing matrix object.
As such, the created vector is not alias safe, and does not take into account
that the underlying matrix memory could be freed (e.g., due to any operation
involving a size change of the matrix).

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> subview1_(const int& n) {
  mat A(n, n, fill::zeros);

  A.submat(0,1,2,3) = randu<mat>(3,3);
  A(span(0,2), span(1,3)) = randu<mat>(3,3);
  A(0,1, size(3,3)) = randu<mat>(3,3);

  mat B = A.submat(0,1,2,3);
  mat C = A(span(0,2), span(1,3) );
  mat D = A(0, 1, size(3,3) );

  A.col(1) = randu<mat>(5,1);
  A(span::all, 1) = randu<mat>(5,1);

  mat X(5, 5, fill::randu);
    
  // get all elements of X that are greater than 0.5
  vec q = X.elem( find(X > 0.5) );
    
  // add 123 to all elements of X greater than 0.5
  X.elem( find(X > 0.5) ) += 123.0;
    
  // set four specific elements of X to 1
  uvec indices = { 2, 3, 6, 8 };
    
  X.elem(indices) = ones<vec>(4);
    
  // add 123 to the last 5 elements of vector a
  vec a(10, fill::randu);
  a.tail(5) += 123.0;
    
  // add 123 to the first 3 elements of column 2 of X
  X.col(2).head(3) += 123;

  return as_doubles_matrix(X);  // Convert from C++ to R
}
```

# Subcube views and slices

A collection of member functions of the `Cube` class that provide subcube views.

## Contiguous views for cube  
      
* `Q.slice(slice_number)`
* `Q.slices(first_slice, last_slice)`
* `Q.row(row_number)`
* `Q.rows(first_row, last_row)`
* `Q.col(col_number)`
* `Q.cols(first_col, last_col)`
* `Q.subcube( first_row, first_col, first_slice, last_row, last_col, last_slice)`
* `Q(span(first_row, last_row), span(first_col, last_col), span(first_slice, last_slice))`
* `Q(first_row, first_col, first_slice, size(n_rows, n_cols, n_slices))`
* `Q(first_row, first_col, first_slice, size(R))` (`R` is a cube)      
* `Q.head_slices(number_of_slices)`
* `Q.tail_slices(number_of_slices)`
* `Q.tube(row, col)`
* `Q.tube(first_row, first_col, last_row, last_col)`
* `Q.tube(span(first_row, last_row), span(first_col, last_col))`
* `Q.tube(first_row, first_col, size(n_rows, n_cols))`

## Non-contiguous views for cube
      
`Q.elem(vector_of_indices)`, `Q(vector_of_indices)`, and
`Q.slices( vector_of_slice_indices)` are instances of `span(a,b)` that can be
replaced by:

* `span()` or `span::all`, to indicate the entire range.
* `span(a)`, to indicate a particular row, column or slice.
  
An individual slice, accessed via `.slice()`, is an instance of the `Mat` class
(a reference to a matrix is provided).

All `.tube()` forms are variants of `.subcube()`, using `first_slice = 0` and
`last_slice = Q.n_slices-1`. The `.tube(row,col)` form uses
`row = first_row = last_row`, and `col = first_col = last_col`.

In the function `Q.elem(vector_of_indices)`, elements specified in
`vector_of_indices` are accessed. `Q` is interpreted as one long vector, with
slice-by-slice and column-by-column ordering of the elements of `Q`. The
`vector_of_indices` must evaluate to a vector of type `uvec` (e.g., generated by
the `find()` function). The aggregate set of the specified elements is treated
as a column vector (e.g., the output of `Q.elem()` is always a column vector).

In the function `Q.slices(vector_of_slice_indices)`, slices specified in
`vector_of_slice_indices` are accessed. The `vector_of_slice_indices` must
evaluate to a vector of type `uvec`.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> subview2_(const int& n) {
  cube A(n, 3, 4, fill::randu);
    
  mat B = A.slice(1); // each slice is a matrix
    
  A.slice(0) = randu<mat>(2,3);
  A.slice(0)(1,2) = 99.0;
    
  A.subcube(0,0,1,  1,1,2)           = randu<cube>(2,2,2);
  A(span(0,1), span(0,1), span(1,2)) = randu<cube>(2,2,2);
  A(0,0,1, size(2,2,2))              = randu<cube>(2,2,2);
    
  // add 123 to all elements of A greater than 0.5
  A.elem( find(A > 0.5) ) += 123.0;
    
  cube C = A.head_slices(2);  // get first two slices
    
  A.head_slices(2) += 123.0;

  mat res = A.slice(0) + B + C.slice(1);

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Subfield views

A collection of member functions of the `field` class that provide subfield
views.      

For a 2D field `F`, the subfields are accessed as:
  
* `F.row(row_number)`
* `F.col(col_number)`
* `F.rows(first_row, last_row)`
* `F.cols(first_col, last_col)`
* `F.subfield(first_row, first_col, last_row, last_col)`
* `F(span(first_row, last_row), span(first_col, last_col))`
* `F(first_row, first_col, size(G))` (`G` is a 2D field)
* `F(first_row, first_col, size(n_rows, n_cols))`

For a 3D field `F`, the subfields are accessed as:
  
* `F.slice(slice_number)`
* `F.slices(first_slice, last_slice)`
* `F.subfield(first_row, first_col, first_slice, last_row, last_col, last_slice)`
* `F(span(first_row, last_row), span(first_col, last_col), span(first_slice, last_slice))`
* `F(first_row, first_col, first_slice, size(G))` (`G` is a 3D field)
* `F(first_row, first_col, first_slice, size(n_rows, n_cols, n_slices))`

Instances of `span(a,b)` can be replaced by:

* `span()` or `span::all`, to indicate the entire range.
* `span(a)`, to indicate a particular row or column.

# Diagonal

`.diag()` is a member functions of `Mat` and `SpMat` with
read/write access to the diagonal in a matrix. The argument can be empty or a
value `k` to specify the diagonal to (`k = 0` by default). The diagonal is
interpreted as a column vector within expressions.

* `k = 0` indicates the main diagonal (default setting)
* `k < 0` indicates the `k`-th sub-diagonal (below main diagonal, towards
  bottom-left corner)
* `k > 0` indicates the `k`-th super-diagonal (above main diagonal, towards
  top-right corner)
      
## Examples

```cpp
[[cpp11::register]] doubles diagonal1_(const int& n) {
  mat X(n, n, fill::randu);

  vec A = X.diag(); // extract the main diagonal
  double B = accu(X.diag(1)); // sum of elements on the first upper diagonal
  double C = accu(X.diag(-1)); // sum of elements on the first lower diagonal

  X.diag() = randu<vec>(n);
  X.diag() += A;
  X.diag() /= B;
  X.diag() *= C;

  sp_mat S = sprandu<sp_mat>(n, n, 0.0);
  S.diag().ones();

  vec v(S.diag());  // copy sparse diagonal to dense vector
  v += X.diag();

  return as_doubles(v);  // Convert from C++ to R
}
```

## Caveat

To calculate only the diagonal elements of a compound expression, use
`diagvec()` or `diagmat()`.

# Each col

`.each_col()` is a member function of `Mat`. It applies a vector operation to
each column of a matrix, and are similar to "broadcasting" in Matlab/Octave.
The argument can be empty, a vector of indices, or a lambda function.
  
| Operation                                | `.each_col()` | `.each_col(vector_of_indices)` | `.each_col(lambda)` |
|------------------------------------------|---------------|--------------------------------|---------------------|
| `+` addition                             | ✓             | ✓                              |                     |
| `+=` in-place addition                   | ✓             | ✓                              |                     |
| `-` subtraction                          | ✓             | ✓                              |                     |
| `-=` in-place subtraction                | ✓             | ✓                              |                     |
| `%` element-wise multiplication          | ✓             | ✓                              |                     |
| `%=` in-place element-wise multiplication| ✓             | ✓                              |                     |
| `/` element-wise division                | ✓             | ✓                              |                     |
| `/=` in-place element-wise division      | ✓             | ✓                              |                     |
| `=` assignment (copy)                    | ✓             | ✓                              |                     |
| `lambda` (lambda function)               |               |                                | ✓                   |
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> each_col1_(const int& n) {
  mat X(n, n + 1, fill::ones);

  // create a vector with n elements ranging from 5 to 10
  vec v = linspace<vec>(5, 10, n);

  // in-place addition of v to each column vector of X
  X.each_col() += v;

  // generate Y by adding v to each column vector of X
  mat Y = X.each_col() + v;

  // subtract v from columns 1 and 2 of X
  X.cols(0, 1).each_col() -= v;

  uvec indices(2);
  indices(0) = 1;
  indices(1) = 2;

  X.each_col(indices) = v;  // copy v to columns 1 and 2 of X

  // lambda function with non-const vector
  X.each_col([](vec& a) { 2 * a; });

  const mat& XX = X;

  // lambda function with const vector
  XX.each_col([](const vec& b) { 3 * b; });

  mat res = X + Y + XX;

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Each row

`.each_row()`, `.each_row(vector_of_indices)`, `.each_row(lambdaction)` are
member functions of `Mat`. These apply a vector operation to each row of
a matrix, and are similar to "broadcasting" in Matlab/Octave.

## Form 1

`.each_row()` supports the following operations:
      
* `+` addition
* `+=` in-place addition
* `-` subtraction
* `-=` in-place subtraction
* `%` element-wise multiplication
* `%=` in-place element-wise multiplication
* `/` element-wise division
* `/=` in-place element-wise division
* `=` assignment (copy)
    
## Form 2

`.each_row(vector_of_indices)` supports the same operations as form 1. The
argument `vector_of_indices` contains a list of indices of the rows to be used,
and it must evaluate to a vector of type `uvec`.
  
## Form 3

`.each_col(lambdaction)` applies the given `lambdaction` to each column
vector. The function must accept a reference to a `Row` object with the same
element type as the underlying matrix.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> each_row1_(const int& n) {
  mat X(n + 1, n, fill::ones);

  // create a vector with n elements ranging from 5 to 10
  rowvec v = linspace<rowvec>(5, 10, n);

  // in-place addition of v to each rows vector of X
  X.each_row() += v;

  // generate Y by adding v to each rows vector of X
  mat Y = X.each_row() + v;

  // subtract v from rows 1 and 2 of X
  X.rows(0, 1).each_row() -= v;

  uvec indices(2);
  indices(0) = 1;
  indices(1) = 2;

  X.each_row(indices) = v;       // copy v to columns 1 and 2 of X

  // lambda function with non-const vector
  X.each_row([](rowvec& a) { a / 2; });

  const mat& XX = X;

  // lambda function with const vector
  XX.each_row([](const rowvec& b) { b / 3; });

  mat res = X + Y + XX;

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Each slice

`.each_slice()` is a member function of `Cube` that applies a matrix operation to
each slice of a cube, with each slice treated as a matrix. It is similar to
"broadcasting" in Matlab/Octave.

## Form 1

`.each_slice(vector_of_indices)`

Supported operations:  
      
* `+` addition
* `+=` in-place addition
* `-` subtraction
* `-=` in-place subtraction 
* `%` element-wise multiplication
* `%=` in-place element-wise multiplication
* `/` element-wise division
* `/=` in-place element-wise division
* `*` matrix multiplication
* `*=` in-place matrix multiplication
* `=` assignment (copy)
   
## Form 2

`.each_slice(lambdaction)`

* The argument _vector_of_indices_ contains a list of indices of the slices to
  be used; it must evaluate to a vector of type `uvec`.
* Arithmetic operations as per form 1 are supported, except for `*` and `*=`
  (e.g., matrix multiplication).

## Form 3

`.each_slice(lambdaction, use_mp)`

* Apply the given `lambdaction` to each slice.
* The function must accept a reference to a `Mat` object with the same element
  type as the underlying cube.

## Form 4

* Apply the given `lambdaction` to each slice, as per form 3.
* The argument `use_mp` is a bool to enable the use of OpenMP for multi-threaded
  execution of `lambdaction` on multiple slices at the same time.
* The order of processing the slices is not deterministic (e.g., slice 2 can be
  processed before slice 1).
* `lambdaction` must be thread-safe, e.g., it must not write to variables
  outside of its scope.
  
## Examples:

```cpp
[[cpp11::register]] doubles_matrix<> each_slice1_(const int& n) {
  cube C(n, n + 1, 6, fill::randu);

  mat M = repmat(linspace<vec>(1, n, n), 1, n + 1);

  C.each_slice() += M;  // in-place addition of M to each slice of C

  cube D = C.each_slice() + M;  // generate D by adding M to each slice of C

  // sum all slices of D into a single n x (n + 1) matrix
  mat D_flat = sum(D, 2);

  uvec indices(2);
  indices(0) = 2;
  indices(1) = 4;

  C.each_slice(indices) = M;  // copy M to slices 2 and 4 in C
  C.each_slice([](mat& X) { X * 2.0; });  // lambda function with non-const matrix
  mat C_flat = sum(C, 2);

  const cube& CC = C;
  CC.each_slice([](const mat& X) { X / 3.0; });  // lambda function with const matrix

  mat CC_flat = sum(CC, 2);

  mat res = C_flat + D_flat + CC_flat;

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Set real

`.set_real(X)` sets the real part of an object. `X` must have the same
size as the recipient object.

## Examples

```cpp
[[cpp11::register]] list set_real1_(const int& n) {
  mat A(n + 1, n - 1, fill::randu);

  cx_mat C(n + 1, n - 1, fill::zeros);

  C.set_real(A);

  return as_complex_matrix(C);  // Convert from C++ to R
}
```
  
## Caveat

To directly construct a complex matrix out of two real matrices, the following
code is faster:

```cpp
[[cpp11::register]] list set_real2_(const int& n) {    
  mat A(n - 1, n + 1, fill::randu);
  mat B(n - 1, n + 1, fill::randu);
  
  cx_mat C = cx_mat(A,B);

  return as_complex_matrix(C);  // Convert from C++ to R
}
```

# Set imaginary

`.set_imaginary(X)` sets the imaginary part of an object. `X` must have the same
size as the recipient object.

## Examples

```cpp
[[cpp11::register]] list set_imag1_(const int& n) {
  mat B(n + 1, n - 1, fill::randu);

  cx_mat C(n + 1, n - 1, fill::zeros);

  C.set_imag(B);

  return as_complex_matrix(C);  // Convert from C++ to R
}
```
  
## Caveat

To directly construct a complex matrix out of two real matrices, the following
code is faster:

```cpp
[[cpp11::register]] list set_imag2_(const int& n) {    
  mat A(n - 1, n + 1, fill::randu);
  mat B(n - 1, n + 1, fill::randu);
  
  cx_mat C = cx_mat(A,B);

  return as_complex_matrix(C);  // Convert from C++ to R
}
```

# Insert columns

`.insert_cols()` is a member function of `Mat`, `Row` and `Cube`. The arguments
can be `colnumber, X` to indicate the column number and the matrix to insert, or
`colnumber, number_of_cols` to indicate the column number and the number of
columns to insert.

The `X` argument inserts a copy of `X` at the specified column. `X` must have
the same number of rows (and slices) as the recipient object.
  
The `number_of_cols` argument expands the object by creating new columns that
are set to zero.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> insert_columns1_(const int& n) {
  mat A(n, n * 2, fill::randu);
  mat B(n, n - 1, fill::ones);

  // at column n - 1, insert a copy of B
  // A will now have 3n - 1 columns
  A.insert_cols(n - 1, B);

  // at column 1, insert 2n zeroed columns
  // B will now have 3n - 1 columns
  B.insert_cols(1, n * 2);

  mat res = A + B;

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Insert rows

`.insert_rows()` is a member function of `Mat`, `Row` and `Cube`. The arguments
can be `rownumber, X` to indicate the row number and the matrix to insert, or
`rownumber, number_of_rows` to indicate the row number and the number of rows to
insert.

The `X` argument inserts a copy of `X` at the specified column. `X` must have
the same number of columns (and slices) as the recipient object.
  
The `number_of_rows` argument expands the object by creating new rows that
are set to zero.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> insert_rows1_(const int& n) {
  mat A(n * 2, n, fill::randu);
  mat B(n - 1, n, fill::ones);

  // at row n - 1, insert a copy of B
  // A will now have 3n - 1 rows
  A.insert_rows(n - 1, B);

  // at row 1, insert 2n zeroed rows
  // B will now have 3n - 1 columns
  B.insert_rows(1, n * 2);

  mat res = A + B;

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Insert slice

`.insert_slices()` is a member function of `Cube`. The arguments can be
`slice_number, X` to indicate the slice number and the matrix to insert, or
`slice_number, number_of_slices` to indicate the slice number and the number of
slices to insert.

The `X` argument inserts a copy of `X` at the specified slice. `X` must have
the same number of columns and rows as the recipient object.
  
The `number_of_slices` argument expands the object by creating new slices that
are set to zero.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> insert_slices1_(const int& n) {
  cube A(n, n, n * 2, fill::randu);
  cube B(n, n, n - 1, fill::ones);

  // At slice n - 1, insert a copy of B
  // A will now have 3n - 1 slices
  A.insert_slices(n - 1, B);

  // At slice 1, insert 2n zeroed slices
  // B will now have 3n - 1 slices
  B.insert_slices(1, n * 2);

  mat res = sum(A + B);

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Shed columns

`.shed_col(row_number)` and `.shed_cols(first_row, last_row)` are member
functions of `Mat`, `Col`, `SpMat`, and `Cube`. With a single scalar argument it
remove the specified column, and with two scalar arguments it removes the
specified range of columns.

`.shed_cols(vector_of_indices)` is a member function of `Mat` and `Col`. With a
vector of indices it must evaluate to a vector of type `uvec` containing the
indices of the columns to remove.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> shed_columns1_(const int& n) {
  mat A(n, n * 5, fill::randu);

  // remove the first column
  A.shed_col(0);

  // remove columns 1 and 2
  A.shed_cols(0, 1);

  // remove columns 2 and 4
  uvec indices(2);
  indices(0) = 1;
  indices(1) = 3;
  A.shed_cols(indices);

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

# Shed rows

`.shed_row(row_number)` and `.shed_rows(first_row, last_row)` are member
functions of `Mat`, `Col`, `SpMat`, and `Cube`. With a single scalar argument it
remove the specified rows, and with two scalar arguments it removes the
specified range of rows.

`.shed_rows(vector_of_indices)` is a member function of `Mat` and `Row`. With a
vector of indices it must evaluate to a vector of type `uvec` containing the
indices of the rows to remove.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> shed_rows1_(const int& n) {
  mat A(n * 5, n, fill::randu);

  // remove the first row
  A.shed_row(0);

  // remove rows 1 and 2
  A.shed_rows(0, 1);

  // remove rows 2 and 4
  uvec indices(2);
  indices(0) = 1;
  indices(1) = 3;
  A.shed_rows(indices);

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

# Shed slices

`.shed_slices()` is a member function of `Cube`. With a single
scalar argument it remove the specified slices, and with two scalar arguments it
removes the specified range of slices. With a vector of indices it must evaluate
to a vector of type `uvec` containing the indices of the rows to remove. The 
arguments can be `slice_number` to indicate the slice number to remove,
`first_slice, last_slice` to indicate the range of slices to remove, or
`vector_of_indices` to indicate the indices of the slices to remove.
  
## Examples

```cpp
[[cpp11::register]] doubles_matrix<> shed_slices1_(const int& n) {
  cube A(n, n, n * 5, fill::randu);

  // remove the first slice
  A.shed_slice(0);

  // remove slices 1 and 2
  A.shed_slices(0, 1);

  // remove slices 2 and 4
  uvec indices(2);
  indices(0) = 1;
  indices(1) = 3;
  A.shed_slices(indices);

  mat res = sum(A, 2);

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

# Swap columns

`.swap_cols( col1, col2 )` is a member functions of `Mat`, `Col`, `Row`, and
`SpMat`. It swaps the contents of the specified columns.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> swap_columns1_(const int& n) {
  mat A(n, n * 5, fill::randu);

  // swap columns 1 and 2
  A.swap_cols(0, 1);

  // swap columns 2 and 4
  A.swap_cols(1, 3);

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

# Swap rows

`.swap_rows( col1, col2 )` is a member functions of `Mat`, `Col`, `Row`, and
`SpMat`. It swaps the contents of the specified rows.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> swap_rows1_(const int& n) {
  mat A(n * 5, n, fill::randu);

  // swap rows 1 and 2
  A.swap_rows(0, 1);

  // swap rows 2 and 4
  A.swap_rows(1, 3);

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

# Swap

`.swap( X )` is a member function of `Mat`, `Col`, `Row`, and `Cube`. It swaps
the contents with object `X`.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> swap1_(const int& n) {
  mat A(n, n + 1, fill::zeros);
  mat B(n * 2, n - 1, fill::ones);

  A.swap(B);

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

# Memory pointer

`.memptr()` is a member function of `Mat`, `Col`, `Row`, and `Cube`. It obtains
a raw pointer to the memory used for storing elements. Data for matrices is
stored in a column-by-column order. Data for cubes is stored in a slice-by-slice
(matrix-by-matrix) order.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> memptr1_(const int& n) {
  mat A(n, n, fill::randu);
  const mat B(n, n, fill::randu);

  double* A_mem = A.memptr();
  const double* B_mem = B.memptr();

  // alter A_mem
  // B_mem is const, so it cannot be altered
  for (int i = 0; i < n * n; ++i) {
    A_mem[i] += 123.0 + B_mem[i];
  }

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

## Caveats

* The pointer becomes invalid after any operation involving a size change or
  aliasing. 
* This function is not recommended for use unless you know what you are doing.

# Column pointer

`.colptr( col_number )` is a member function of the `Mat` class that obtains a
raw pointer to the memory used by elements in the specified column.

## Examples

```cpp
[[cpp11::register]] doubles_matrix<> colptr1_(const int& n) {
  mat A(n, n, fill::randu);

  // pointer to the memory of the first column of A
  double* Acol1_mem = A.colptr(0);

  // alter memory
  for (int i = 0; i < n; ++i) {
    Acol1_mem[i] += 123.0;
  }

  return as_doubles_matrix(A);  // Convert from C++ to R
}
```

## Caveats

* The pointer becomes invalid after any operation involving a size change or
  aliasing.
* This function is not recommended for use unless you know what you are doing.
* It is safer to use `submat()` instead.

# Iterators

Iterators for traverse over all elements within the specified range. These
return the column/row/slice of an object as a `uword` type.

## Member functions

Dense matrices and vectors (`Mat`, `Col`, and `Row`):

* `.begin()` is an iterator referring to the first element.
* `.end()` is an iterator referring to the past the end element.
* `.begin_col(col_number)` is an iterator referring to the first element of the
  specified column.
* `.end_col(col_number)` is an iterator referring to the past-the-end element of
  the specified column.
* `begin_row(row_number)` is an iterator referring to the first element of the
  specified row.
* `end_row(row_number)` is an iterator referring to the past-the-end element of
  the specified row.

Cubes (`Cube`): 
      
* `begin()` is an iterator referring to the first element.
* `end()` is an iterator referring to the past-the-end element.
* `begin_slice(slice_number)` iterator referring to the first element of the
  specified slice.
* `end_slice(slice_number)` iterator referring to the past-the-end element of
  the specified slice.

Sparse matrices (`SpMat`):
    
* `begin()` is an iterator referring to the first element.
* `end()` is an iterator referring to the past-the-end element.
* `begin_col(col_number)` is an iterator referring to the first element of the
  specified column.
* `end_col(col_number)` is an iterator referring to the past-the-end element of
  the specified column.
* `begin_row(row_number)` is an iterator referring to the first element of the
  specified row.
* `end_row(row_number)` is an iterator referring to the past-the-end element of
  the specified row.

Dense submatrices and subcubes (`submatrix` and `subcube`):

* `span(row, col)` and `span(row, col, slice)` can be used to specify the range
  of elements to iterate over.

## Iterator types

Dense matrices and vectors (`Mat`, `Col`, and `Row`):
    
* `mat::iterator`, `vec::iterator` and `rowvec::iterator` are random access
  iterators, for read/write access to elements (which are stored column by
  column).
* `mat::const_iterator`, `vec::const_iterator` and `rowvec::const_iterator`
  are random access iterators, for read-only access to elements (which are
  stored column by column)
* `mat::col_iterator`, `vec::col_iterator` and `rowvec::col_iterator` random
  access iterators, for read/write access to the elements of specified columns.
* `mat::const_col_iterator`, `vec::const_col_iterator` and
  `rowvec::const_col_iterator` are random access iterators, for read-only access
  to the elements of specified columns.
* `mat::row_iterator` is a bidirectional iterator, for read/write access to the
  elements of specified rows.
* `mat::const_row_iterator` is a bidirectional iterator, for read-only access to
  the elements of specified rows.
* `vec::row_iterator` and `rowvec::row_iterator` are random access iterators,
  for read/write access to the elements of specified rows.
* `vec::const_row_iterator` and `rowvec::const_row_iterator` are random access
  iterators, for read-only access to the elements of specified rows.

Cubes (`Cube`):

* `cube::iterator` is a random access iterator, for read/write access to
  elements. The elements are ordered slice by slice; the elements within each
  slice are ordered column by column.
* `cube::const_iterator` is a random access iterator, for read-only access to
  elements.
* `cube::slice_iterator` is a random access iterator, for read/write access to
  the elements of a particular slice. The elements are ordered column by column.
* `cube::const_slice_iterator` is a random access iterator, for read-only access
  to the elements of a particular slice.

Sparse matrices (`SpMat`):

* `sp_mat::iterator` is a bidirectional iterator, for read/write access to
  elements (which are stored column by column).
* `sp_mat::const_iterator` is a bidirectional iterator, for read-only access to
  elements (which are stored column by column).
* `sp_mat::col_iterator` is a bidirectional iterator, for read/write access to
  the elements of a specific column.
* `sp_mat::const_col_iterator` is a bidirectional iterator, for read-only access
  to the elements of a specific column.
* `sp_mat::row_iterator` is a bidirectional iterator, for read/write access to
  the elements of a specific row.
* `sp_mat::const_row_iterator` is a bidirectional iterator, for read-only access
  to the elements of a specific row.

## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> iterators1_(const int& n) {
  mat X(n, n + 1, fill::randu);

  mat::iterator it = X.begin();
  mat::iterator it_end = X.end();

  for (; it != it_end; ++it) {
    (*it) += 123.0;
  }

  mat::col_iterator col_it = X.begin_col(1);    // start of column 1
  mat::col_iterator col_it_end = X.end_col(n);  //   end of column n

  for (; col_it != col_it_end; ++col_it) {
    (*col_it) = 321.0;
  }

  return as_doubles_matrix(X);  // Convert from C++ to R
}
```

```cpp
[[cpp11::register]] doubles_matrix<> iterators2_(const int& n) {
  cube X(n, n + 1, n + 2, fill::randu);

  cube::iterator it = X.begin();
  cube::iterator it_end = X.end();

  for (; it != it_end; ++it) {
    (*it) += 123.0;
  }

  cube::slice_iterator s_it = X.begin_slice(1);    // start of slice 1
  cube::slice_iterator s_it_end = X.end_slice(n);  // end of slice n

  for (; s_it != s_it_end; ++s_it) {
    (*s_it) = 321.0;
  }

  mat res = sum(X, 2);

  return as_doubles_matrix(res);  // Convert from C++ to R
}
```

```cpp
[[cpp11::register]] doubles_matrix<> iterators3_(const int& n) {
  sp_mat X = sprandu<sp_mat>(n, n * 2, 0.1);

  sp_mat::iterator it = X.begin();
  sp_mat::iterator it_end = X.end();

  for (; it != it_end; ++it) {
    (*it) += 123.0;
  }

  return as_doubles_matrix(X);  // Convert from C++ to R
}
```

```cpp
[[cpp11::register]] doubles_matrix<> iterators4_(const int& n) {
  mat X(n, n, fill::randu);

  for (double& val : X(span(0, 1), span(1, 1))) {
    val = 123.0;
  }

  return as_doubles_matrix(X);  // Convert from C++ to R
}
```

## Caveats

* Writing a zero value into a sparse matrix through an iterator will
  invalidate all current iterators associated with the sparse matrix.
* To modify the non-zero elements in a safer manner, use `.transform()` or
  `.for_each()` instead of iterators.
* For `submatrix` and `subcube` the iterators are intended only to be used with
  range-based for loops. Any other use is not supported. For example, the direct
  use of the `.begin()` and `.end()` functions, as well as the underlying
  iterators types is not supported. The implementation of submatrices and
  subcubes uses short-lived temporary objects that are subject to automatic
  deletion, and as such are error-prone to handle manually.

# Compatibility container functions

Member functions for the `Col` and `Row` classes to mimic the functionality of
containers in the C++ standard library:  
      
* `.front()` accesses the first element in a vector.
* `.back()` accesses the last element in a vector.

Member functions for the `Col`, `Row`, `Mat`, `Cube` and `SpMat` classes to
mimic the functionality of containers in the C++ standard library:

* `.clear()` removes the elements from an object.
* `.empty()` returns `true` if the object has no elements and `false` if the
  object has one or more elements.    
* `.size()` returns the total number of elements in an object.

## Examples

```cpp
[[cpp11::register]] doubles compatibility1_(const int& n) {
  vec X(n, fill::randu);

  writable::doubles res = {X.front(), X.back()};

  res.attr("names") = strings({"front", "back"});

  return res;
}
```

```cpp
[[cpp11::register]] integers compatibility2_(const int& n) {
  mat X(n, n, fill::randu);

  writable::integers res(2);
  res[0] = X.n_rows;

  X.clear();
  res[1] = X.n_rows;

  res.attr("names") = strings({"before", "after"});

  return res;
}
```

# Convert matrix to column

`.as_col()` is a member function of the `Mat` class, it returns a flattened
version of the matrix as a column vector. Flattening is done by concatenating
all columns.
  
## Examples

```cpp
[[cpp11::register]] doubles as_col1_(const int& n) {
  mat M(n, n + 1, fill::randu);
  vec V = M.as_col();
  return as_doubles(V);
}
```

# Convert matrix to row

`.as_row()` is a member function of the `Mat` class, it returns a flattened
version of the matrix as a row vector. Flattening is done by concatenating
all rows.
  
## Examples

```cpp
[[cpp11::register]] doubles as_row1_(const int& n) {
  mat M(n, n + 1, fill::randu);
  vec V = M.as_row();
  return as_doubles(V);
}
```

## Caveat

Converting columns to rows is faster than converting rows to columns.

# Convert column to matrix

`.col_as_mat(col_number)` is a member function of the `Cube` class, it returns a
matrix of the specified cube column and the number of rows is preserved. Given
a cube of size `R x C x S`, the resultant matrix size is `R x S`.

## Examples

```cpp
[[cpp11::register]] list col_as_mat1_(const int& n) {
  cube C(n, n + 1, n + 2, fill::randu);
  mat M = C.col_as_mat(0);  // size n x (n + 1)
  
  writable::list res(5);
  res[0] = as_doubles_matrix(C.slice(0));
  res[1] = as_doubles_matrix(C.slice(1));
  res[2] = as_doubles_matrix(C.slice(2));
  res[3] = as_doubles_matrix(C.slice(3));
  res[4] = as_doubles_matrix(M);

  res.attr("names") = strings({"slice0", "slice1", "slice2", "slice3",
    "col_as_mat"});

  return res;
}
```

# Convert column to matrix

`.row_as_mat(row_number)` is a member function of the `Cube` class, it returns a
matrix of the specified cube row and the number of columns is preserved. Given
a cube of size `R x C x S`, the resultant matrix size is `S x C`.

## Examples

```cpp
[[cpp11::register]] list row_as_mat1_(const int& n) {
  cube C(n, n + 1, n + 2, fill::randu);
  mat M = C.row_as_mat(0);  // size (n + 2) x (n + 1)

  writable::list res(5);
  res[0] = as_doubles_matrix(C.slice(0));
  res[1] = as_doubles_matrix(C.slice(1));
  res[2] = as_doubles_matrix(C.slice(2));
  res[3] = as_doubles_matrix(C.slice(3));
  res[4] = as_doubles_matrix(M);

  res.attr("names") = strings({"slice0", "slice1", "slice2", "slice3",
    "row_as_mat"});

  return res;
}
```

# Convert sparse matrix to dense matrix

`.as_dense()` is a member function of the `SpMat` class, it avoids the
construction of an intermediate sparse matrix representation of the expression.

## Examples
    
```cpp
[[cpp11::register]] doubles as_dense1_(const int& n) {
  sp_mat A;
  A.sprandu(n, n, 0.1);

  // extract column 1 of A directly into dense column vector
  colvec c = A.col(0).as_dense();

  // store the sum of each column of A directly in dense row vector
  rowvec r = sum(A).as_dense();

  return as_doubles(c + r.t());
}
```

# Dense matrix and vector transposition

`.t()` is a member function of the `Mat`, `Col` and `Row` classes, it returns a
transposed copy of the object. For real matrices, the transpose is a simple
transposition of the elements. For complex matrices, the transpose is a
Hermitian conjugate transposition of the elements (e.g., the signs of the
imaginary components are flipped).
  
## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> transpose1_(const int& n) {
  mat A(n, n + 1, fill::randu);
  mat B = A.t();
  return as_doubles_matrix(B);
}
```

# Sparse matrix transposition

`.st()` is a member function of the `SpMat` classe, it returns a
transposed copy of the object. For real matrices, it is not applicable. For
complex matrices, the transpose is a simple transposition of the elements (e.g.,
the signs of imaginary components are not flipped).
  
## Examples
    
```cpp
[[cpp11::register]] doubles_matrix<> transpose2_(const int& n) {
  sp_mat A;
  A.sprandu(n, n + 1, 0.1);
  sp_mat B = A.t();
  return as_doubles_matrix(B);
}
```

# Matrix inversion

`.i()` is a member function of the `Mat` class, it provides an inverse of the
matrix. If the matrix is not square sized, a `std::logic_error` exception is
thrown. If the matrix appears to be singular, the output matrix is reset and a
`std::runtime_error` exception is thrown.
  
##  Examples
    
```cpp
[[cpp11::register]] doubles inverse1_(const doubles_matrix<>& a,
                                      const doubles b) {
  mat A = as_Mat(a);
  vec B = as_Col(b);

  mat X = inv(A);
  vec Y = X * B;

  return as_doubles(Y);
}
```

## Caveats

* If the matrix is known to be symmetric positive definite, `inv_sympd()`.
* To solve a system of linear equations, such as `Z = inv(X) * Y`, `solve()` can
  be faster and/or more accurate.

# Maximum and minimum

`.min()` and `.max()` are member functions of the `Mat`, `Col`, `Row`, and
`Cube` classes. These return the minimum and maximum values of the object,
respectively. For objects with complex numbers, absolute values are used for
comparison.
  
## Examples
    
```cpp
[[cpp11::register]] doubles maxmin1_(const int& n) {
  mat A = randu<mat>(n, n);

  writable::doubles res(2);
  res[0] = A.max();
  res[1] = A.min();

  res.attr("names") = strings({"max", "min"});

  return res;
}
```

# Linear index of maximum and minimum

`.index_min()` and `.index_max()` are member functions of the `Mat`, `Col`,
`Row`, and `Cube` classes. They return the linear index of the minimum and
maximum values of the object, respectively. For objects with complex numbers,
absolute values are used for comparison. The returned index is of type `uword`.

## Examples
    
```cpp
[[cpp11::register]] doubles index_maxmin1_(const int& n) {
  mat A = randu<mat>(n, n);

  writable::doubles res(6);
  res[0] = static_cast<int>(A.index_max());
  res[1] = static_cast<int>(A.index_min());
  res[2] = A(0, 0);
  res[3] = A(1, 0);
  res[4] = A(0, 1);
  res[5] = A(1, 1);

  res.attr("names") = strings({"index_max", "index_min", "element0", "element1",
    "element2", "element3"});

  return res;
}
```

# In-range

`.in_range(** i **)` is a member function of `Mat`, `Col`, `Row`, `Cube`,
`SpMat` and `field`, it returns `true` if the given location or span is
currently valid and `false` if the object is empty, the location is out of
bounds, or the span is out of bounds.

| Function                                                                                      | Mat | Col | Row | Cube | SpMat | Field |
|-----------------------------------------------------------------------------------------------|-----|-----|-----|------|-------|-------|
| `.in_range(span(start, end))`                                                                 | ✓   | ✓   | ✓   | ✓    | ✓     | ✓     |
| `.in_range(row, col)`                                                                         | ✓   | ✓   | ✓   | ✓    | ✓     | ✓     |
| `.in_range(span(start_row, end_row), span(start_col, end_col))`                               | ✓   | ✓   | ✓   |      | ✓     | ✓     |
| `.in_range(row, col, slice)`                                                                  |     |     |     | ✓    |       | ✓     |
| `.in_range(span(start_row, end_row), span(start_col, end_col), span(start_slice, end_slice))` |     |     |     | ✓    |       | ✓     |
| `.in_range(first_row, first_col, size(X))` (`X` is a matrix or field)                         | ✓   | ✓   | ✓   |      | ✓     | ✓     |
| `.in_range(first_row, first_col, size(n_rows, n_cols))`                                       | ✓   | ✓   | ✓   |      | ✓     | ✓     |
| `.in_range(first_row, first_col, first_slice, size(Q))` (`Q` is a cube or field)              |     |     |     | ✓    |       | ✓     |
| `.in_range(first_row, first_col, first_slice, size(n_rows, n_cols, n_slices))`                |     |     |     | ✓    |       | ✓     |

Instances of `span(a,b)` can be replaced by:

* `span()` or `span::all` to indicate the entire range.
* `span(a)` to indicate a particular row, column, or slice.
  
## Examples
    
```cpp
[[cpp11::register]] logicals in_range1_(const int& n) {
  mat A(n, n + 1, fill::randu);

  writable::logicals res(3);
  res[0] = A.in_range(0, 0);
  res[1] = A.in_range(3, 4);
  res[2] = A.in_range(4, 5);

  res.attr("names") = strings({"in_range00", "in_range34", "in_range45"});

  return res;
}
```

# Is empty

`.is_empty()` is a member function of the `Mat`, `Col`, `Row`, `Cube`, `SpMat`,
and `field` classes. It returns `true` if the object has no elements and `false`
if the object has one or more elements.

## Examples

```cpp
[[cpp11::register]] logicals is_empty1_(const int& n) {
  mat A(n, n + 1, fill::randu);

  writable::logicals res(2);
  res[0] = A.is_empty();

  A.reset();
  res[1] = A.is_empty();

  res.attr("names") = strings({"before_reset", "after_reset"});

  return res;
}
```

# Is vector/column vector/row vector

`.is_vec()`, `.is_colvec()` and `.is_rowvec()` are member functions of `Mat`
and `SpMat`.

* `.is_vec()` returns `true` if the matrix can be interpreted as a vector
  (either column or row vector) and `false` otherwise.
* `.is_colvec()` returns `true` if the matrix can be interpreted as a column
  vector and `false` otherwise.  
* `.is_rowvec()` returns `true` if the matrix can be interpreted as a row
  vector and `false` otherwise.
    
## Examples

```cpp
[[cpp11::register]] logicals is_vec1_(const int& n) {
  mat A(n, 1, fill::randu);
  mat B(1, n, fill::randu);
  mat C(0, 1, fill::randu);
  mat D(1, 0, fill::randu);

  writable::logicals res(5);
  res[0] = A.is_vec();
  res[1] = A.is_colvec();
  res[2] = B.is_rowvec();
  res[3] = C.is_colvec();
  res[4] = D.is_rowvec();

  res.attr("names") = strings({"Nx1_is_vec", "Nx1_is_colvec", "1xN_is_rowvec",
    "0x1_is_colvec", "1x0_is_rowvec"});

  return res;
}
```

## Caveat

Do not assume that the vector has elements if these functions return `true`. It
is possible to have an empty vector (e.g., 0x1 as in the examples).

# Is sorted

`.is_sorted()`, `.is_sorted(sort_direction)` and
`.is_sorted(sort_direction, dim)` are member function of `Mat`, `Row`, and
`Col`. For matrices and vectors with complex numbers, order is checked via
absolute values.

If the object is a vector, these return a `bool` indicating whether the elements
are sorted. If the object is a matrix, these return a `bool` indicating whether
the elements are sorted in each column (`dim = 0`, default) or each row
(`dim = 1`), and the `dim` argument is optional.

The `sort_direction` argument is optional, `sort_direction` can be one of the
following strings:
    
* `"ascend"`: the elements are ascending, consecutive elements can be equal, and
  this is the default operation.
* `"descend"`: the elements are descending, and consecutive elements can be
  equal.
* `"strictascend"`: the elements are strictly ascending, and consecutive
  elements cannot be equal.
* `"strictdescend"`: the elements are strictly descending, and consecutive
  elements cannot be equal.

## Examples

```cpp
[[cpp11::register]] logicals is_sorted1_(const int& n) {
  vec a(n, fill::randu);
  vec b = sort(a);
  mat A(10, 10, fill::randu);

  writable::logicals res(4);
  res[0] = a.is_sorted();
  res[1] = b.is_sorted();
  res[2] = A.is_sorted("descend", 1);
  res[4] = A.is_sorted("ascend", 1);

  res.attr("names") = strings({"a_sorted", "b_sorted", "A_descend",
    "A_ascend"});

  return res;
}
```

# Is upper triangular/lower triangular

`.is_trimatu()` and `.is_trimatl()` are member functions of `Mat` and `SpMat`.
`.is_trimatu()` returns `true` if the matrix is upper triangular (e.g., the
matrix is square sized and all elements below the main diagonal are zero) and
`false` otherwise. `.is_trimatl()` returns `true` if the matrix is lower
triangular (e.g., the matrix is square sized and all elements above the main
diagonal are zero) and `false` otherwise.
  
## Examples

```cpp
[[cpp11::register]] logicals is_triangular1_(const int& n) {
  mat A(n, n, fill::randu);
  mat B = trimatl(A);

  writable::logicals res(3);
  res[0] = B.is_trimatu();
  res[1] = B.is_trimatl();

  B.reset();
  res[2] = B.is_trimatu();

  res.attr("names") = strings({"is_trimatu", "is_trimatl",
    "is_trimatu_after_reset"});

  return res;
}
```

## Caveat

If these functions return `true`, do not assume that the matrix contains
non-zero elements on or above/below the main diagonal. It is possible to have an
empty matrix (e.g., 0x0 as in the examples).

# Is diagonal

`is_diagmat()` is a member function of `Mat` and `SpMat`. It returns `true` if
the matrix is diagonal (e.g., all elements outside of the main diagonal are
zero). If the matrix is not square sized, a `std::logic_error` exception is
thrown.

## Examples

```cpp
[[cpp11::register]] logicals is_diagonal1_(const int& n) {
  mat A(n, n, fill::randu);
  mat B = diagmat(A);

  writable::logicals res(3);
  res[0] = A.is_diagmat();
  res[1] = B.is_diagmat();

  A.reset();
  res[2] = A.is_diagmat();

  res.attr("names") = strings({"A_diagmat", "B_diagmat",
    "A_diagmat_after_reset"});

  return res;
}
```

## Caveat

If this function returns `true`, do not assume that the matrix contains non-zero
elements on the main diagonal only. It is possible to have an empty matrix
(e.g., 0x0 as in the examples).

# Is square

`.is_square()` is a member function of the `Mat` and `SpMat` classes. It returns
`true` if the matrix is square sized (e.g., the number of rows is equal to the
number of columns) and `false` otherwise.

## Examples

```cpp
[[cpp11::register]] logicals is_square1_(const int& n) {
  mat A(n, n, fill::randu);
  mat B = diagmat(A);

  writable::logicals res(3);
  res[0] = A.is_square();
  res[1] = B.is_square();

  A.reset();
  res[2] = A.is_square();

  res.attr("names") = strings({"A_square", "B_square",
    "A_square_after_reset"});

  return res;
}
```

## Caveats

If this function returns `true`, do not assume that the matrix is non-empty. It
is possible to have an empty matrix (e.g., 0x0 as in the examples).

# Is symmetric

`.is_symmetric()` is a member function of the `Mat` and `SpMat` classes. It
returns `true` if the matrix is symmetric (e.g., the matrix is square sized and
the transpose is equal to the original matrix) and `false` otherwise.

## Examples

```cpp
[[cpp11::register]] logicals is_symmetric1_(const int& n) {
  mat A(n, n, fill::randu);
  mat B = symmatu(A);

  writable::logicals res(3);
  res[0] = A.is_symmetric();
  res[1] = B.is_symmetric();

  A.reset();
  res[2] = A.is_symmetric();

  res.attr("names") = strings({"A_symmetric", "B_symmetric",
    "A_symmetric_after_reset"});

  return res;
}
```

## Caveats

If this function returns `true`, do not assume that the matrix is non-empty. It
is possible to have an empty matrix (e.g., 0x0 as in the examples).

# Is hermitian

`.is_hermitian()` is a member function of the `Mat` and `SpMat` classes. It
returns `true` if the matrix is Hermitian or self-adjoint (e.g., the matrix is
square sized and the conjugate transpose is equal to the original matrix) and
`false` otherwise.

## Examples

```cpp
[[cpp11::register]] logicals is_hermitian1_(const int& n) {
  cx_mat A(n, n, fill::randu);
  cx_mat B = A.t() * A;

  writable::logicals res(3);
  res[0] = A.is_hermitian();
  res[1] = B.is_hermitian();

  A.reset();
  res[2] = A.is_hermitian();

  res.attr("names") = strings({"A_hermitian", "B_hermitian",
    "A_hermitian_after_reset"});

  return res;
}
```

## Caveats

If this function returns `true`, do not assume that the matrix is non-empty. It
is possible to have an empty matrix (e.g., 0x0 as in the examples).

# Is symmetric/hermitian positive definite

`.is_sympd()` and `.is_sympd(tol)` are a member function of the `Mat` and
`SpMat` classes. It returns `true` if the matrix is symmetric/hermitian
positive definite within a tolerance (e.g., the matrix is square sized and all
its eigenvalues are positive) and `false` otherwise. The `tol` argument is
optional, the default is `tol = 100 * datum::eps * norm(X, "fro")`.

## Examples

```cpp
[[cpp11::register]] logicals is_sympd1_(const int& n) {
  mat A(n, n, fill::randu);
  mat B = A * A.t();

  writable::logicals res(3);
  res[0] = A.is_sympd();
  res[1] = B.is_sympd();

  A.reset();
  res[2] = A.is_sympd();

  res.attr("names") = strings({"A_sympd", "B_sympd",
    "A_sympd_after_reset"});

  return res;
}
```

# Is zero

`.is_zero()` and `.is_zero(tol)` are a member function of the `Mat`, `Col`,
`Row`, `Cube`, and `SpMat` classes. It returns `true` if all elements are zero
within a tolerance and `false` otherwise. For complex numbers, each component
(real and imaginary) is checked separately. The `tol` argument is optional.

## Examples

```cpp
[[cpp11::register]] logicals is_zero1_(const int& n) {
  mat A(n, n, fill::randu);
  cube B(n, n, n, fill::zeros);
  sp_mat C(n, n);

  writable::logicals res(3);
  res[0] = A.is_zero(0.005);
  res[1] = B.is_zero(0.005);
  res[2] = C.is_zero(0.005);

  res.attr("names") = strings({"A_is_zero", "B_is_zero", "C_is_zero"});

  return res;
}
```

# Is finite

`.is_finite()` is a member function of the `Mat`, `Col`, `Row`, `Cube`, and
`SpMat` classes. It returns `true` if all elements are finite and `false`
otherwise.

## Examples

```cpp
[[cpp11::register]] logicals is_finite1_(const int& n) {
  mat A(n, n, fill::randu);
  cube B(n, n, n, fill::randu);
  sp_mat C(n, n);

  // Insert infinite values
  B(0, 0, 0) = datum::inf;
  C(0, 0) = -1.0 * datum::inf;

  writable::logicals res(3);
  res[0] = A.is_finite();
  res[1] = B.is_finite();
  res[2] = C.is_finite();

  res.attr("names") = strings({"A_is_finite", "B_is_finite", "C_is_finite"});

  return res;
}
```

# Has infinity

`.has_inf()` is a member function of the `Mat`, `Col`, `Row`, `Cube`, and
`SpMat` classes. It returns `true` if the object contains at least one infinite
value and `false` otherwise.

## Examples

```cpp
[[cpp11::register]] logicals has_inf1_(const int& n) {
  mat A(n, n, fill::randu);
  cube B(n, n, n, fill::randu);
  sp_mat C(n, n);

  // Insert infinite values
  B(0, 0, 0) = datum::inf;
  C(0, 0) = -1.0 * datum::inf;

  writable::logicals res(3);
  res[0] = A.has_inf();
  res[1] = B.has_inf();
  res[2] = C.has_inf();

  res.attr("names") = strings({"A_has_inf", "B_has_inf", "C_has_inf"});

  return res;
}
```

# Has not-a-number

`.has_nan()` is a member function of the `Mat`, `Col`, `Row`, `Cube`, and
`SpMat` classes. It returns `true` if the object contains at least one
not-a-number (NaN) value and `false` otherwise.

## Examples

```cpp
[[cpp11::register]] logicals has_nan1_(const int& n) {
  mat A(n, n, fill::randu);
  cube B(n, n, n, fill::randu);
  sp_mat C(n, n);

  // Insert NaN values
  B(0, 0, 0) = datum::nan;
  C(0, 0) = -1.0 * datum::nan;

  writable::logicals res(3);
  res[0] = A.has_nan();
  res[1] = B.has_nan();
  res[2] = C.has_nan();

  res.attr("names") = strings({"A_has_nan", "B_has_nan", "C_has_nan"});

  return res;
}
```

## Caveat

`NaN` is not equal to anything, even itself.
