API — ICoreBlocks/ICoreMath/Foundation
The public contract of 12 header(s) under src/ICoreBlocks/ICoreMath/Foundation — 12 class/struct definition(s), 424 declaration(s). Each section shows the header's banner and its public (and protected-virtual) surface exactly as the file writes it.
ICoreComplexMatrix.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreComplexMatrix.h
ICoreComplexMatrix#
ICoreComplexMatrix.h:32 · class · pImpl · 37 declaration(s)
A 2-D matrix of complex doubles -- the console's Kind::Complex value (the core-MATLAB row C0.2), and the wrapper that replaces the N x 2 [real, imag] convention every complex-producing function...
class ICoreComplexMatrix {
public:
ICoreComplexMatrix();
ICoreComplexMatrix(size_t rows, size_t columns);
explicit ICoreComplexMatrix(const std::complex<double>& scalar);
ICoreComplexMatrix(const ICoreComplexMatrix& other);
ICoreComplexMatrix& operator=(const ICoreComplexMatrix& other);
~ICoreComplexMatrix();
// Every real matrix IS a complex one with zero imaginary parts; the two
// builders below are the only way in, so a half-built value cannot exist.
static ICoreComplexMatrix fromReal(const ICoreMatrix& real);
static ICoreComplexMatrix fromParts(const ICoreMatrix& real, const ICoreMatrix& imaginary,
std::string* whyNot = nullptr);
static ICoreComplexMatrix fromList(const std::vector<ICoreComplexVariable>& values);
size_t size1() const;
size_t size2() const;
size_t numel() const;
bool isScalar() const;
bool isEmpty() const;
bool isVector() const;
std::complex<double>& operator()(size_t i, size_t j);
const std::complex<double>& operator()(size_t i, size_t j) const;
// The four projections MATLAB names (C5.44-C5.46), each answering a REAL
// matrix of the same shape.
ICoreMatrix realPart() const;
ICoreMatrix imagPart() const;
ICoreMatrix magnitude() const; // abs(z)
ICoreMatrix phase() const; // angle(z), atan2(imag, real)
// MATLAB's isreal() is a question about the STORAGE, not the values: a
// value built as complex answers false even when every imaginary part is
// zero. hasImaginaryPart() is the other question -- whether any imaginary
// part is actually non-zero -- and it is what decides whether a result
// narrows back to a real matrix.
bool hasImaginaryPart() const;
ICoreMatrix narrowToReal() const; // the real parts, whatever the imaginary ones are
ICoreComplexMatrix conjugate() const;
ICoreComplexMatrix transpose() const; // .' -- no conjugation
ICoreComplexMatrix ctranspose() const; // ' -- conjugate transpose
ICoreComplexMatrix negate() const;
// Element-wise arithmetic, broadcasting a 1x1 operand the way ICoreMatrix
// does. A size mismatch fills `whyNot` and answers a default matrix.
ICoreComplexMatrix add(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
ICoreComplexMatrix subtract(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
ICoreComplexMatrix elementWiseMultiply(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
ICoreComplexMatrix elementWiseDivide(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
ICoreComplexMatrix elementWisePower(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
ICoreComplexMatrix multiply(const ICoreComplexMatrix& other, std::string* whyNot = nullptr) const;
// A \ b with complex entries: a square system solved exactly, a tall one
// in the least-squares sense -- the same two readings ICoreMatrix::solve
// gives, over a column-pivoting QR.
ICoreComplexMatrix leftDivide(const ICoreComplexMatrix& rhs, std::string* whyNot = nullptr) const;
ICoreComplexMatrix elementWiseSqrt() const;
ICoreComplexMatrix elementWiseExp() const;
ICoreComplexMatrix elementWiseLog() const;
double normFrobenius() const;
bool isApprox(const ICoreComplexMatrix& other, const double& tolerance) const;
// The PRINCIPAL matrix square root and the general real matrix power
// (the core-MATLAB row C7.17), both over the complex field because that
// is the only field they are total in: [1 2; 3 4] has a negative
// eigenvalue, so its square root is complex and a real-only routine can
// only answer NaN -- which is what this console did until this row.
//
// Both go through Eigen's Schur-based matrix functions rather than an
// eigendecomposition, and the difference is not academic: a DEFECTIVE
// matrix has too few eigenvectors to diagonalise, so an eig-based root
// divides by a singular V and answers noise. [1 1; 0 1] is the smallest
// one, its square root is exactly [1 0.5; 0 1], and MATLAB answers that
// (measured against R2026a, 2026-09-02) because MATLAB's sqrtm is
// Schur-based too.
//
// `matrixPower` is MATLAB's A^p for a non-integer p, and p = 0.5 agrees
// with matrixSquareRoot to the last few bits on both sides. A singular
// matrix with a negative exponent has no answer and says so.
ICoreComplexMatrix matrixSquareRoot(std::string* whyNot = nullptr) const;
ICoreComplexMatrix matrixPower(const double& exponent, std::string* whyNot = nullptr) const;
// ---- the rows that PRODUCE a complex value -----------------------------
// MATLAB's fft/ifft along `dim` (1 = down the columns, 2 = across the
// rows, 0 = the first non-singleton), padded or truncated to `n` when n is
// non-zero. C10.1, C10.2; fft2/ifft2 (C10.3) are the two applied in turn,
// which is what MATLAB's own fft2 is.
static ICoreComplexMatrix dft(const ICoreComplexMatrix& x, size_t n, int dim,
const bool& inverse, std::string* whyNot = nullptr);
// Eigenvalues as a COLUMN, and MATLAB's [V, D] and [V, D, W] forms
// (C7.8). `vectors`/`leftVectors` are filled only when asked for.
static bool eigen(const ICoreMatrix& a,
ICoreComplexMatrix& values,
ICoreComplexMatrix* vectors,
ICoreComplexMatrix* leftVectors,
std::string* whyNot = nullptr);
// The generalized problem A*v = lambda*B*v (C7.8's two-argument form).
static bool eigenGeneralized(const ICoreMatrix& a, const ICoreMatrix& b,
ICoreComplexMatrix& values,
ICoreComplexMatrix* vectors,
std::string* whyNot = nullptr);
private:
class Impl; // the two-line residue; state lives here
std::unique_ptr<Impl> impl;
};ICoreComplexVariable.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreComplexVariable.h
////////////////////////// / Scalar operations /////////////////////////
ICoreComplexVariable#
ICoreComplexVariable.h:7 · class · pImpl · 32 declaration(s)
class ICoreComplexVariable {
public:
ICoreComplexVariable();
ICoreComplexVariable(const double& real, const double& imaginary);
static ICoreComplexVariable fromPolar(const double& magnitude, const double& angle);
////////////////////////////
/// Scalar operations
///////////////////////////
ICoreComplexVariable operator*(double scalar) const;
ICoreComplexVariable operator/(double scalar) const;
////////////////////////////
/// Basic arithmetic operators
///////////////////////////
ICoreComplexVariable operator+(const ICoreComplexVariable &other) const;
ICoreComplexVariable operator-(const ICoreComplexVariable &other) const;
ICoreComplexVariable operator*(const ICoreComplexVariable &other) const;
ICoreComplexVariable operator/(const ICoreComplexVariable &other) const;
////////////////////////////
/// Conjugate
///////////////////////////
ICoreComplexVariable conjugate() const;
////////////////////////////
/// Other operations
///////////////////////////
double magnitude() const;
double magnitudeSquared() const;
double phase() const;
double complexDistance(const ICoreComplexVariable &other) const;
ICoreComplexVariable normalized() const;
ICoreComplexVariable exp() const;
ICoreComplexVariable log() const;
ICoreComplexVariable pow(const double& exponent) const;
bool isApprox(const ICoreComplexVariable &other, const double& tol = 1e-12) const;
bool isZero(double tol) const;
bool isReal(double tol) const;
ICoreComplexVariable rotate(const double &angle) const;
////////////////////////////
/// Setters/Getters
///////////////////////////
void setReal(const double& newReal);
void setImaginary(const double& newImaginary);
double getReal() const;
double getImaginary() const;
std::string getVariableAsString() const;
friend std::ostream& operator<<(std::ostream& os, const ICoreComplexVariable& mat);
// Complex values live in std::vector (eigenvalues, root loci), so the type
// stays copyable; the unique_ptr residue below deletes the implicit copy
// operations, so all four are written out in the .cpp.
ICoreComplexVariable(const ICoreComplexVariable& other);
ICoreComplexVariable& operator=(const ICoreComplexVariable& other);
ICoreComplexVariable(ICoreComplexVariable&& other) noexcept;
ICoreComplexVariable& operator=(ICoreComplexVariable&& other) noexcept;
~ICoreComplexVariable();
private:
class Impl; // the two-line residue; state lives here
std::unique_ptr<Impl> impl;
};ICoreFactorUpdates.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreFactorUpdates.h
ICoreFactorUpdates#
ICoreFactorUpdates.h:35 · class · 1 declaration(s)
Core MATLAB's matfun/ factorization UPDATES and the two norm estimators (the core-MATLAB board, C7.30): given a factorization you already have, the factorization of a slightly different matrix, c...
class ICoreFactorUpdates {
public:
// ---- planerot ------------------------------------------------------
// MATLAB's own Givens rotation, and it is here in the header's contract
// because the four QR updates below are DEFINED by which rotations they
// apply in which order: G = [x'; -x(2) x(1)]/r with r = norm(x), and the
// IDENTITY (not a sign-flipped rotation) when x(2) is already zero. That
// last case is why a zero row does not change the factors at all.
static void planeRotation(double a, double b, double& c, double& s, double& r);
// ---- cholupdate ----------------------------------------------------
// R is the UPPER Cholesky factor of A (MATLAB's `chol`), and this answers
// the upper factor of A + x*x' (update) or A - x*x' (downdate).
//
// LINPACK's `dchud`/`dchdd` rather than a transcription: MATLAB's
// cholupdate is a built-in with no `.m` to read, so this is the published
// algorithm those built-ins are named for and NOT a line-for-line copy of
// what MATLAB runs. It agrees with MATLAB to the factorization band, and
// that is the claim -- not the bit-for-bit agreement the transcribed rows
// on this board make.
//
// `failedRank` is MATLAB's second output `p`: 0 when the answer is a full
// Cholesky factor, and otherwise the 1-based step at which A - x*x' stopped
// being positive definite. A downdate that fails leaves `out` untouched.
static bool choleskyUpdate(const ICoreMatrix& r, const ICoreMatrix& x, bool downdate,
ICoreMatrix& out, size_t* failedRank = nullptr,
std::string* whyNot = nullptr);
// ---- qrupdate ------------------------------------------------------
// A = Q*R and this answers the factors of A + u*v'. Same footing as
// cholupdate: the published rank-1 QR update, not a transcription.
static bool qrUpdate(const ICoreMatrix& q, const ICoreMatrix& r,
const ICoreMatrix& u, const ICoreMatrix& v,
ICoreMatrix& qOut, ICoreMatrix& rOut,
std::string* whyNot = nullptr);
// ---- qrinsert / qrdelete -------------------------------------------
// The factors of A with one column or row inserted before index `j`
// (1-based, MATLAB's own indexing), or deleted at it. `byRow` picks
// MATLAB's 'row' orientation over its default 'col'.
//
// These two ARE transcriptions -- `qrinsert.m` and `qrdelete.m` carry a
// readable path beside their built-in call, and it is the rotations in
// that path, in that order, that are written out here.
static bool qrInsert(const ICoreMatrix& q, const ICoreMatrix& r,
size_t j, const ICoreMatrix& x, bool byRow,
ICoreMatrix& qOut, ICoreMatrix& rOut,
std::string* whyNot = nullptr);
static bool qrDelete(const ICoreMatrix& q, const ICoreMatrix& r,
size_t j, bool byRow,
ICoreMatrix& qOut, ICoreMatrix& rOut,
std::string* whyNot = nullptr);
// ---- normest -------------------------------------------------------
// `normest.m` transcribed: the 2-norm by power iteration to a RELATIVE
// tolerance, default 1e-6. It is an ESTIMATE and not a value -- MATLAB's
// own answer differs from `norm(A, 2)` in the sixth or seventh digit by
// design -- so a parity case over it is a BAND and never a comparison, and
// the console's row says so. `iterations` is MATLAB's second output.
//
// The first iterate is the vector of COLUMN 1-norms, which is what makes
// the iteration reproducible: nothing here is random unless a product
// comes out exactly zero, and MATLAB reaches for `rand` only there.
static bool norm2Estimate(const ICoreMatrix& a, double tol,
double& estimate, size_t* iterations = nullptr,
std::string* whyNot = nullptr);
// ---- cdf2rdf / rsf2csf ---------------------------------------------
// The two conversions between the eigen-decomposition's REAL and COMPLEX
// spellings, and they run in opposite directions:
//
// cdf2rdf takes `[V, D] = eig(A)` for a real A -- complex conjugate
// eigenvalues down D's diagonal -- and answers a REAL pair with those
// eigenvalues in 2 x 2 blocks. `cdf2rdf.m` transcribed.
//
// rsf2csf takes `[U, T] = schur(A)`, the REAL Schur form with 2 x 2
// blocks on T's diagonal, and answers the COMPLEX Schur form with the
// eigenvalues on the diagonal. `rsf2csf.m` transcribed.
//
// cdf2rdf assumes conjugate pairs are ADJACENT, which is what eig answers,
// and says so rather than producing nonsense when they are not.
static bool complexToRealBlockDiagonal(const ICoreComplexMatrix& v,
const ICoreComplexMatrix& d,
ICoreMatrix& vOut, ICoreMatrix& dOut,
std::string* whyNot = nullptr);
static bool realToComplexSchur(const ICoreMatrix& u, const ICoreMatrix& t,
ICoreComplexMatrix& uOut, ICoreComplexMatrix& tOut,
std::string* whyNot = nullptr);
// ---- condest -------------------------------------------------------
// An estimate of `cond(A, 1)` -- `norm(A, 1) * norm(inv(A), 1)` with the
// second factor ESTIMATED rather than computed, so no inverse is formed.
//
// **This is Hager's algorithm, and MATLAB's condest is not.** MATLAB calls
// `normest1`, the Higham-Tisseur BLOCK generalisation of Hager, whose
// default block size is 2 and whose starting matrix has a first column of
// 1/n and the rest **random** (`normest1.m`). So MATLAB's own condest is
// not reproducible run to run on a matrix large enough for the estimate to
// matter, which is exactly what C0.9 keeps off this console. The single-
// vector Hager iteration here is deterministic and finds the EXACT 1-norm
// on small matrices -- which is where the two agree, and the row says so
// rather than claiming agreement it cannot have at size.
//
// `witness` is MATLAB's second output `v`: a vector with `norm(A*v, 1)`
// small relative to `norm(v, 1)`, normalised so `norm(v, 1) == 1`.
static bool conditionEstimate(const ICoreMatrix& a, double& estimate,
ICoreMatrix* witness = nullptr,
std::string* whyNot = nullptr);
// ---- polyeig -------------------------------------------------------
// The polynomial eigenvalue problem
// (A0 + lambda*A1 + ... + lambda^p * Ap) * x = 0
// by MATLAB's own linearization (`polyeig.m`): build the n*p by n*p pair
// (A, B) whose generalized eigenvalues are the answer, and hand it to the
// QZ algorithm C7.8 already wired up. `vectors`, when asked for, is filled
// with MATLAB's own choice -- for each eigenvalue, the block of the big
// eigenvector with the smallest normalized residual, scaled to unit
// 2-norm.
static bool polynomialEigenvalues(const std::vector<ICoreMatrix>& coefficients,
ICoreComplexMatrix& values,
ICoreComplexMatrix* vectors,
std::string* whyNot = nullptr);
};
};ICoreGeneralizedSVD.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreGeneralizedSVD.h
ICoreGeneralizedSVD#
ICoreGeneralizedSVD.h:43 · class · 0 declaration(s)
MATLAB's gsvd -- the generalized singular value decomposition of a pair of matrices with the same number of columns (the core-MATLAB board, C7.30's last name, and its own module because it is the...
class ICoreGeneralizedSVD {
public:
// [U, V, X, C, S] = gsvd(A, B), or gsvd(A, B, "econ") for `economy`,
// which needs m >= p or n >= p and answers U and V with at most p
// columns (MATLAB's own restriction and shape). B must have A's column
// count; anything else says so and answers false.
static bool decompose(const ICoreMatrix& a, const ICoreMatrix& b, bool economy,
ICoreMatrix& u, ICoreMatrix& v, ICoreMatrix& x,
ICoreMatrix& c, ICoreMatrix& s,
std::string* whyNot = nullptr);
// sigma = gsvd(A, B): the generalized singular values as a column, which
// is sqrt(diag(C'*C) ./ diag(S'*S)) read off the same decomposition --
// MATLAB's one-output form, computed by the same path so the two
// spellings cannot disagree.
// MATLAB accepts the "econ" word here too and reads the values off the
// economy decomposition; `economy` is that.
static bool values(const ICoreMatrix& a, const ICoreMatrix& b, bool economy,
ICoreMatrix& sigma, std::string* whyNot = nullptr);
};
};ICoreLogic.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreLogic.h
ICoreLogic#
ICoreLogic.h:32 · class · 2 declaration(s)
Core MATLAB's yes/no answers (the core-MATLAB board, C9.1-C9.3, C9.7) and the three non-finite tests that answer the same shape (C5.48).
class ICoreLogic {
public:
// isnan / isinf / isfinite, element-wise (C5.48).
enum class FiniteTest { IsNaN, IsInf, IsFinite };
static ICoreMatrix finiteTest(const FiniteTest& which, const ICoreMatrix& a,
std::string* whyNot = nullptr);
// any(A), any(A, dim), any(A, "all") -- and all(...) in the same forms
// (C9.1). `dim` is 1 (down columns), 2 (across rows) or 0 for MATLAB's
// default, the first non-singleton dimension; `everyElement` is the "all"
// word, which reduces the whole array to one answer.
static ICoreMatrix anyOf(const ICoreMatrix& a, int dim, bool everyElement,
std::string* whyNot = nullptr);
static ICoreMatrix allOf(const ICoreMatrix& a, int dim, bool everyElement,
std::string* whyNot = nullptr);
// The six comparisons, element-wise, with MATLAB's rules (C1.5): a scalar
// operand broadcasts, any other size mismatch is refused the way MATLAB
// refuses it, and `==` is EXACT -- the console's eq() compared within
// 1e-12 until this landed, so `0.1 + 0.2 == 0.3` answered 1 where MATLAB
// answers 0 (C0.1 chose MATLAB's meaning for a core name). A NaN is
// unequal to everything including itself, which falls out of IEEE and is
// the one place these differ from the connectives below: `<` and `==`
// take a NaN, `and`/`or` refuse it.
enum class Comparison { Eq, Ne, Lt, Le, Gt, Ge };
static ICoreMatrix compare(const Comparison& which, const ICoreMatrix& a,
const ICoreMatrix& b, std::string* whyNot = nullptr);
// and / or / xor, and not (C9.2). A scalar operand broadcasts; any other
// size mismatch is refused, as MATLAB refuses it. `a & b`, `a | b` and
// `~a` are the same three operations under MATLAB's operator spelling
// since C1.6 retired `|`/`&` as concatenation, and they reach exactly
// this code -- the function form and the operator cannot drift.
enum class Connective { And, Or, Xor };
static ICoreMatrix connective(const Connective& which, const ICoreMatrix& a,
const ICoreMatrix& b, std::string* whyNot = nullptr);
static ICoreMatrix negate(const ICoreMatrix& a, std::string* whyNot = nullptr);
// isequal(A, B, ...) and isequaln (C9.3): every value must have the same
// shape AND the same entries. They differ on one point only -- isequal
// says NaN is not equal to NaN (it follows ==), isequaln says it is.
static bool equalValues(const std::vector<ICoreMatrix>& values, bool nanEqualsNan,
std::string* whyNot = nullptr);
// issorted(v), issorted(A), issorted(A, "rows") (C9.7). Without "rows" a
// matrix is sorted when each COLUMN is non-decreasing and a row vector when
// it is non-decreasing along its length (the first non-singleton rule
// again); with it, the rows must be non-decreasing lexicographically.
static bool isSorted(const ICoreMatrix& a, bool byRows, std::string* whyNot = nullptr);
};
};ICoreMatrix.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreMatrix.h
Two members below (toEigen2DMatrixXd / fromEigen) name Eigen::MatrixXd in their signatures, so this header cannot be parsed without Eigen. It used to arrive for free from the force-included SDK-root pch.h; that Eigen block now lives in ICoreMath's own umbrella and is named here explicitly.
ICoreMatrix#
ICoreMatrix.h:49 · class · pImpl · 263 declaration(s)
class ICoreMatrix {
public:
////////////////////////////////////////////////////
///
/// Constructors and Base Methods
///
////////////////////////////////////////////////////
ICoreMatrix();
ICoreMatrix(size_t d1, size_t d2, size_t d3 = 1);
explicit ICoreMatrix(const double& scalar);
explicit ICoreMatrix(const std::vector<double>& vec1d);
explicit ICoreMatrix(const std::vector<std::vector<double>>& vec2d);
explicit ICoreMatrix(const std::vector<std::vector<std::vector<double>>>& vec3d);
static ICoreMatrix identity(const size_t& n);
static ICoreMatrix zeros(const size_t &size1, const size_t &size2);
static ICoreMatrix zerosVector(const size_t &size1);
static ICoreMatrix ones(const size_t &size1, const size_t &size2);
static ICoreMatrix onesVector(const size_t &size1);
// Both are ROW vectors, MATLAB's shape (the core-MATLAB board, C4.9 and
// C4.10); they were columns until 2026-09-02. n = 1 is the END point, and
// logspace sweeps to pi itself, not to 10^pi, when the upper exponent is
// exactly pi -- MATLAB's documented special case.
static ICoreMatrix linspace(const double& a, const double& b, const size_t& n);
static ICoreMatrix logspace(const double& a, const double& b, const size_t& n);
static ICoreMatrix rand(const size_t& size1, const size_t& size2); // uniform [0,1)
static ICoreMatrix randn(const size_t& size1, const size_t& size2); // standard normal
// This is the setter method
double& operator()(size_t i, size_t j, size_t k = 0);
// This is the getter method
double operator()(size_t i, size_t j, size_t k = 0) const;
size_t size1() const;
size_t size2() const;
size_t size3() const;
std::vector<double> rawData() const;
void assignRawDataDirectly(const size_t& size1, const size_t& size2, const size_t& size3, const std::vector<double>& newRawDataVector);
size_t numberOfRows() const;
size_t numberOfColumns() const;
void resize(size_t d1, size_t d2 = 1, size_t d3 = 1);
void overwriteIndex(const double& newVal, const size_t i, const size_t j = 0, const size_t k = 0);
void cloneValueFromMatrix(const ICoreMatrix* other);
void cloneValueFromMatrix(const ICoreMatrix& other);
void clear();
////////////////////////////////////////////////////
///
/// Verifications
///
////////////////////////////////////////////////////
bool verifySize(const size_t &d1, const size_t &d2 = 1, const size_t &d3 = 1) const;
bool verifySquareMatrix() const;
bool verifyElementWiseInteger() const;
bool verifyElementWiseNotNegative() const;
bool verifyElementWiseStrictlyPositive() const;
bool verifyZerosMatrix(double tol = 1e-12) const;
bool verifySalarZeroMatrix(double tol) const;
bool verifyOnesMatrix(double tol = 1e-12) const;
bool verifyValidVectorSize() const;
bool isScalar() const;
// Shape and structure predicates (the core-MATLAB board, C4.18 and C7.24).
// Each answers about the LANGUAGE-level shape, which for this console is
// always 2-D: isMatrix() is therefore true for everything it can hold, and
// isEmpty() is false for everything it can hold -- both are answers, not
// stubs, and both stop being trivial the day C3.9 gives the console an
// empty matrix.
// MATLAB's LOGICAL class (decision C0.3, rows C3.6 and C9.4). A flag rather
// than a second value kind, because that is the shape of MATLAB's own rule:
// a logical array holds ones and zeros like any other, and what differs is
// the CLASS it reports and what the class-reading operations do with it.
//
// isLogical() is this matrix a logical array
// asLogical() the same numbers, said to be (or not to be) logical --
// for the operations that ANSWER one (a comparison, a
// connective, any/all, every is* question) and for the few
// that carry one through unchanged
// logicalOf() MATLAB's logical(A): non-zero becomes 1, and a NaN is
// REFUSED rather than converted ("NaN values cannot be
// converted to logicals"), which is MATLAB's own sentence
//
// Arithmetic needs no rule here: every operator builds a fresh matrix and
// a fresh matrix is not logical, which is exactly MATLAB's behaviour
// (`true + 1` is the double 2).
bool isLogical() const;
ICoreMatrix asLogical(bool yes = true) const;
static ICoreMatrix logicalOf(const ICoreMatrix& a, std::string* whyNot = nullptr);
bool isEmpty() const;
bool isVector() const; // 1 x n or n x 1
bool isRow() const; // 1 x n
bool isColumn() const; // n x 1
bool isMatrix() const; // 2-D, i.e. anything this class holds today
bool isDiagonal() const; // every off-diagonal entry exactly 0
bool isUpperTriangular() const;
bool isLowerTriangular() const;
// Every non-zero entry within `lower` bands below and `upper` above the
// main diagonal. isDiagonal() is isBanded(0, 0).
bool isBanded(const size_t& lower, const size_t& upper) const;
bool isSymmetric(double tol = 1e-9) const;
bool isPositiveDefinite() const;
////////////////////////////////////////////////////
///
/// Conversions to std vectors
///
////////////////////////////////////////////////////
// Convert full 3D matrix to nested vectors
std::vector<std::vector<std::vector<double>>> toVector3D() const;
std::vector<std::vector<double>> toVector2D() const;
// Get a 2D slice by fixing dim1 index (like a "matrix" at i)
std::vector<std::vector<double>> getSlice2D(const size_t& i) const;
// // Get a 1D column vector fixing dim1 and dim3 (returns all j for fixed i,k)
// std::vector<double> getColumn(const size_t& i, const size_t& k = 0) const;
// // Get a 1D row vector fixing dim1 and dim2 (returns all k for fixed i,j)
// std::vector<double> getRow(const size_t& i, const size_t& j) const;
ICoreMatrix getColumn(const size_t& j) const;
ICoreMatrix getRow(const size_t& i) const;
// Get a 1D vector of the entire entries for the entire matrix
std::vector<double>& getRawDataVector();
std::vector<double> getRawDataVector_Copy() const;
std::pair<bool, ICoreMatrix> sliceMatrix2D(const size_t &row_start, const size_t &row_end, const size_t &col_start,
const size_t &col_end) const;
////////////////////////////////////////////////////
///
/// Algebra
///
////////////////////////////////////////////////////
// To Apply operators
void operator+=(const ICoreMatrix& other);
ICoreMatrix operator+(const ICoreMatrix& other) const;
void operator-=(const ICoreMatrix& other);
ICoreMatrix operator-(const ICoreMatrix& other) const;
ICoreMatrix operator*(const double& scale) const;
friend ICoreMatrix operator*(double scale, const ICoreMatrix& m);
ICoreMatrix operator*(const ICoreMatrix& other) const;
ICoreMatrix operator|(const ICoreMatrix& other) const;
ICoreMatrix operator&(const ICoreMatrix& other) const;
// Element-wise addition: this + other
ICoreMatrix& add(const ICoreMatrix& other);
// Element-wise subtraction: this - other
ICoreMatrix& subtract(const ICoreMatrix& other);
// Product
double dot(const ICoreMatrix& other) const;
ICoreMatrix cross(const ICoreMatrix& other) const;
ICoreMatrix multiply(const ICoreMatrix& other) const;
ICoreMatrix& multiplyElementWiseByScalar(const double& factor);
ICoreMatrix multiplyElementWiseByMatrix(ICoreMatrix& other) const;
ICoreMatrix transpose() const;
size_t rank(const double& tol = 1e-12) const;
// Inverse
ICoreMatrix inverse() const;
ICoreMatrix syntraInverse() const;
// Eigen-backed scalar quantities
double determinant() const; // square 2D only
double trace() const; // square 2D only
double conditionNumber() const; // ratio of largest/smallest singular value
// Eigen-backed factorizations / generalized inverses
ICoreMatrix pseudoInverse(double tol = -1.0) const; // Moore-Penrose (SVD); tol<0 => auto
ICoreMatrix sqrtm() const; // principal matrix square root (square)
ICoreMatrix choleskyL() const; // lower factor L, A = L*L^T (square SPD)
ICoreMatrix nullSpace() const; // basis for ker(A) as columns
// LDLT decomposition (symmetric, possibly indefinite): P^T*L*diag(D)*L^T*P = A
ICoreMatrix ldltL() const;
ICoreMatrix ldltD() const;
ICoreMatrix ldltP() const;
// Matrix exponential
ICoreMatrix expm() const;
ICoreMatrix logm() const;
// Norms
double norm_Frobenius() const;
double norm_RMS() const;
double norm_1() const;
// The largest singular value, from the SVD. It was a power iteration with
// an iteration count and a 1e-6 tolerance until 2026-09-02: MATLAB's
// norm(A, 2) is exact, so a parity case over it could only ever be run in
// the iterative band (the core-MATLAB board, C7.5).
double norm_2() const;
double norm_inf() const;
double norm_nuclear() const; // sum of singular values, 2D only
ICoreMatrix concatHorizontal(const ICoreMatrix& other) const;
ICoreMatrix concatVertical(const ICoreMatrix& other) const;
std::vector<ICoreComplexVariable> eigenvalues() const;
ICoreMatrix eigenvectors() const; // real part; exact for symmetric matrices
// The SELF-ADJOINT eigendecomposition, and the only pair of accessors here
// whose two halves are in one order: entry k of symmetricEigenvalues() is
// the eigenvalue of column k of symmetricEigenvectors(), both ASCENDING,
// which is the order MATLAB's `eig` answers a symmetric matrix in.
//
// ⚠ eigenvalues() AND eigenvectors() CANNOT BE PAIRED. The first runs
// Eigen's general EigenSolver and the second, for a symmetric matrix, the
// self-adjoint one; the two order their answers by different rules, so
// `eigenvalues()[k]` is not in general the eigenvalue of column k of
// `eigenvectors()`. Anything that needs the pairing -- cmdscale's B, a
// covariance's principal axes -- asks for it here.
ICoreMatrix symmetricEigenvalues() const; // n x 1, ascending
ICoreMatrix symmetricEigenvectors() const; // n x n, column k for value k
// A*x = rhs. Square: the exact solution. Tall: the least-squares one
// (column-pivoting QR, MATLAB's backslash). Wide: refused.
ICoreMatrix solve(const ICoreMatrix &rhs) const;
// QR decomposition (thin): A = Q*R
ICoreMatrix qrQ() const;
ICoreMatrix qrR() const;
// LU decomposition with partial pivoting (square only): P*A = L*U
ICoreMatrix luL() const;
ICoreMatrix luU() const;
ICoreMatrix luP() const;
// SVD (thin): A = U*diag(S)*V^T; svdS is a column vector of singular values
ICoreMatrix svdU() const;
ICoreMatrix svdS() const;
ICoreMatrix svdV() const;
// ---- MATLAB's factorizations, in the shapes MATLAB answers (C7.9-C7.13) ----
//
// The six thin factors above are the console's own extras (qrq, svdu, ...)
// and keep their meaning. MATLAB's `[Q, R] = qr(A)` is the FULL
// factorization -- Q is m x m and R is m x n -- so these are separate
// methods rather than a change of shape under the old names, and MATLAB's
// economy form (`qr(A, 0)`, `svd(A, 'econ')`) is the flag.
ICoreMatrix qrFactorQ(const bool& economy = false) const;
ICoreMatrix qrFactorR(const bool& economy = false) const;
// The COLUMN-PIVOTED QR, A*P = Q*R: MATLAB's three-output form. A
// different factorization, not a rearrangement of the one above.
ICoreMatrix qrPivotedQ() const;
ICoreMatrix qrPivotedR() const;
ICoreMatrix qrPivotedP() const;
// A = U*S*V', with S the m x n MATRIX MATLAB answers rather than the
// vector svdS() gives.
ICoreMatrix svdFactorU(const bool& economy = false) const;
ICoreMatrix svdFactorS(const bool& economy = false) const;
ICoreMatrix svdFactorV(const bool& economy = false) const;
// The UPPER Cholesky factor, A = R'*R -- MATLAB's default, where
// choleskyL() is the lower one. `failedPivot` answers MATLAB's `p`: 0 when
// A is positive definite, otherwise the 1-based pivot at which the
// factorization failed, and the result is then the leading (p-1) x (p-1)
// block that did succeed. Hand-written rather than Eigen's LLT for exactly
// that reason: a decomposition that only says "it failed" cannot answer p.
ICoreMatrix choleskyUpper(size_t* failedPivot = nullptr) const;
// MATLAB's two null-space bases. nullSpaceOrthonormal() is null(A): the
// columns of V whose singular values fall at or below
// max(size(A)) * eps(largest) -- an ORTHONORMAL basis, where nullSpace()
// above answers Eigen's kernel, which is not orthonormal. Both answer a
// matrix with ZERO columns when A has full column rank.
// nullSpaceRational() is null(A, 'r'): read off the reduced row echelon
// form, so its entries are the small rationals a hand calculation gives.
ICoreMatrix nullSpaceOrthonormal() const;
ICoreMatrix nullSpaceRational() const;
// ---- the three similarity transforms MATLAB names (C7.21-C7.23) --------
//
// Each answers BOTH halves of its transform, because the one-output form
// is not checkable without the other half:
// hess: A = P*H*P' P orthogonal, H upper Hessenberg
// schur: A = U*T*U' U orthogonal, T the REAL Schur form (quasi-
// upper-triangular: a 1x1 block per real
// eigenvalue, a 2x2 per complex pair)
// balance: B = T\A*T T a permutation times a diagonal of exact
// powers of two, so no rounding is introduced
// A non-square argument reports through ICoreRunDiagnosis and answers the
// default matrix, the way every factorization here does.
ICoreMatrix hessenbergH() const;
ICoreMatrix hessenbergP() const;
ICoreMatrix schurT() const;
ICoreMatrix schurU() const;
// balanced() is MATLAB's B and balancingSimilarity() its T. `permute` is
// MATLAB's default; false is balance(A, "noperm"), which scales without
// reordering.
ICoreMatrix balanced(const bool& permute = true) const;
ICoreMatrix balancingSimilarity(const bool& permute = true) const;
// diag(): vector -> diagonal matrix; 2D matrix -> its diagonal as a column vector
ICoreMatrix diag() const;
// Kronecker product, 2D only
ICoreMatrix kron(const ICoreMatrix& other) const;
// Element-wise math (operates over all entries, including 3D)
ICoreMatrix elementWiseAbs() const;
ICoreMatrix elementWiseExp() const;
ICoreMatrix elementWiseLog() const;
ICoreMatrix elementWiseSqrt() const;
ICoreMatrix elementWiseSin() const;
ICoreMatrix elementWiseCos() const;
ICoreMatrix elementWiseTan() const;
ICoreMatrix elementWisePow(const double& exponent) const;
ICoreMatrix elementWiseMin(const ICoreMatrix& other) const; // broadcasts a scalar operand
ICoreMatrix elementWiseMax(const ICoreMatrix& other) const; // broadcasts a scalar operand
ICoreMatrix clamp(const double& lo, const double& hi) const;
// Core MATLAB's elementary element-wise math, over all entries including
// 3D. The binary forms broadcast a scalar
// operand exactly as elementWiseMin/Max do; a size mismatch reports
// through ICoreRunDiagnosis and returns the default matrix.
ICoreMatrix elementWiseSign() const;
ICoreMatrix elementWiseFloor() const;
ICoreMatrix elementWiseCeil() const;
ICoreMatrix elementWiseFix() const; // round toward zero
ICoreMatrix elementWiseRound(const int& digits = 0) const; // half AWAY from zero, not banker's
ICoreMatrix elementWiseRoundSignificant(const int& digits) const;
ICoreMatrix elementWiseSec() const;
ICoreMatrix elementWiseCsc() const;
ICoreMatrix elementWiseCot() const;
ICoreMatrix elementWiseSinh() const;
ICoreMatrix elementWiseCosh() const;
ICoreMatrix elementWiseTanh() const;
ICoreMatrix elementWiseDeg2Rad() const;
ICoreMatrix elementWiseSpacing() const; // eps(x): distance to the next double
ICoreMatrix elementWiseRad2Deg() const;
ICoreMatrix elementWiseHypot(const ICoreMatrix& other) const;
ICoreMatrix elementWiseAtan2(const ICoreMatrix& other) const; // atan2(this, other) = atan2(y, x)
ICoreMatrix elementWiseMod(const ICoreMatrix& other) const; // sign follows the DIVISOR
ICoreMatrix elementWiseRem(const ICoreMatrix& other) const; // sign follows the DIVIDEND
// Integer-domain element-wise math. Each reports a non-integer (or, for
// factorial, a negative) entry through whyNot and returns the default
// matrix, the way the Sylvester family reports a singular operator --
// MATLAB raises an error on the same inputs, so returning a plausible
// number would be the one answer that cannot be right.
ICoreMatrix elementWiseFactorial(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseGcd(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseLcm(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
// The prime trio of the core-MATLAB board (C5.31), which belongs with the
// integer-domain family above because it refuses on the same grounds.
// primes(n) is a ROW vector of every prime up to n, by the sieve;
// isprime tests each entry (and MATLAB accepts only non-negative whole
// numbers, so -3 and 2.5 are refused rather than answered false);
// primeFactors(n) is MATLAB's factor(), ascending with repeats, a SCALAR
// argument only. factor(1) is 1 and factor(0) is 0, both as MATLAB has it.
static ICoreMatrix primes(const double& upTo, std::string* whyNot = nullptr);
ICoreMatrix elementWiseIsPrime(std::string* whyNot = nullptr) const;
static ICoreMatrix primeFactors(const double& n, std::string* whyNot = nullptr);
// divisors(n) -- every POSITIVE divisor of n, ascending, as a row.
// MATLAB's is a Symbolic Math Toolbox name that answers a double for a
// double argument, which is why it sits here beside the prime trio rather
// than only on the symbolic side (the toolbox board's T4.18).
// divisors(0) is 0 and the sign of n is ignored, both measured.
static ICoreMatrix divisors(const double& n, std::string* whyNot = nullptr);
// The real-only family: MATLAB leaves the reals for these inputs and the
// real* names exist precisely to refuse instead. nthroot is real by
// definition and needs an odd integer root for a negative base.
ICoreMatrix elementWiseRealSqrt(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseRealLog(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseRealPow(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseNthRoot(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
// expm1/log1p keep the digits exp(x)-1 and log(1+x) lose near zero.
// log1p refuses below -1 for the reason the real* family exists.
ICoreMatrix elementWiseExpm1() const;
ICoreMatrix elementWiseLog1p(std::string* whyNot = nullptr) const;
// log2/log10, the two remaining bases (the core-MATLAB board, C5.7, C5.8).
// Total maps like elementWiseLog(): a negative entry is NaN here, and the
// caller that means MATLAB's log2/log10 asks hasNegativeEntry() first --
// MATLAB answers a COMPLEX number there and this console has no complex
// kind yet (C0.2), so sqrt/log/log2/log10 refuse rather than answer NaN.
ICoreMatrix elementWiseLog2() const;
ICoreMatrix elementWiseLog10() const;
bool hasNegativeEntry() const; // a negative zero is NOT negative, as IEEE has it
// [f, e] = log2(x), MATLAB's two-output form (C5.7) -- and it is std::frexp
// entry by entry: f in [0.5, 1) with x = f * 2^e. It is TOTAL on the reals
// where one-output log2 is not, negatives included (log2(-6) is complex,
// but [f, e] = log2(-6) is -0.75 and 3), so neither answer refuses.
// Zero, the infinities and NaN all report e = 0 and f = x, MATLAB's own.
ICoreMatrix elementWiseLog2Mantissa() const;
ICoreMatrix elementWiseLog2Exponent() const;
// meshgrid/ndgrid (C4.15). One builder for both, because ndgrid IS
// meshgrid transposed: meshgrid(x, y) is numel(y) BY numel(x) with x
// running along the row, ndgrid(x, y) is numel(x) by numel(y) with x
// running down the column. `second` picks the output (X/Y, or the ndgrid
// pair); the orientation of x and y themselves is ignored, as MATLAB
// ignores it.
static ICoreMatrix grid(const ICoreMatrix& x, const ICoreMatrix& y,
const bool& second, const bool& transposed,
std::string* whyNot = nullptr);
// Trigonometry in DEGREES (the core-MATLAB board, C5.21). These are not
// sin(deg2rad(x)): the angle never becomes radians at a multiple of 90, so
// sind(180) is exactly +0, cosd(90) is exactly +0 and tand(90) is Inf --
// which is what makes secd/cscd/cotd answer Inf there rather than 1e16.
// The inverse three refuse outside [-1, 1], where MATLAB goes complex.
ICoreMatrix elementWiseSinDeg() const;
ICoreMatrix elementWiseCosDeg() const;
ICoreMatrix elementWiseTanDeg() const;
ICoreMatrix elementWiseSecDeg() const;
ICoreMatrix elementWiseCscDeg() const;
ICoreMatrix elementWiseCotDeg() const;
ICoreMatrix elementWiseAsinDeg(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAcosDeg(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAtanDeg() const;
// The same three, one number at a time. The Aerospace conversions
// (the toolbox board's T7.10 and T7.11) work in degrees and are
// written scalar-at-a-time over a caller's rows, and MATLAB's own
// geodesy calls sind/cosd/atan2d rather than a radian route -- so the
// exactness above is part of THEIR answers too, and a second copy of the
// octant reduction beside them would be a copy that eventually disagrees.
// atan2d is the pair's inverse and is exact on the axes: atan2Degrees(1, 0)
// is exactly 90, which atan2(1, 0)*180/pi is not.
static double sinDegrees(const double& degrees);
static double cosDegrees(const double& degrees);
static double atan2Degrees(const double& y, const double& x);
// The inverse trigonometry and the whole hyperbolic family (C5.15, C5.17,
// C5.19, C5.20). Each reciprocal is 1/f and each inverse-reciprocal is
// f(1/x) -- MATLAB's own composition, measured bitwise, and the one that
// carries the signed zero out to the pole: acot(-0) is -pi/2 and coth(-0)
// is -Inf. The seven that take a whyNot refuse where MATLAB goes complex,
// with the boundary itself inside the domain: atanh(1), acoth(1) and
// asech(0) are all Inf, and asech(-0) is NaN.
ICoreMatrix elementWiseAsin(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAcos(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAtan() const;
ICoreMatrix elementWiseAsec(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAcsc(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAcot() const;
ICoreMatrix elementWiseSech() const;
ICoreMatrix elementWiseCsch() const;
ICoreMatrix elementWiseCoth() const;
ICoreMatrix elementWiseAsinh() const;
ICoreMatrix elementWiseAcosh(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAtanh(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAsech(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseAcsch() const;
ICoreMatrix elementWiseAcoth(std::string* whyNot = nullptr) const;
// The error function and its inverses (C5.40). erfinv/erfcinv answer NaN
// outside their domain rather than refusing, because MATLAB does.
ICoreMatrix elementWiseErf() const;
ICoreMatrix elementWiseErfc() const;
ICoreMatrix elementWiseErfInv() const;
ICoreMatrix elementWiseErfcInv() const;
// The gamma and beta families (C5.41, C5.42). gamma() is Inf at every
// non-positive integer; gammaln/beta/betaln refuse a negative entry and
// betainc refuses an x outside [0, 1], each exactly where MATLAB errors.
// The `upper` flag is MATLAB's 'upper' tail option.
ICoreMatrix elementWiseGamma() const;
ICoreMatrix elementWiseGammaLn(std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseGammaInc(const ICoreMatrix& a, const bool& upper,
std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseBeta(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseBetaLn(const ICoreMatrix& other, std::string* whyNot = nullptr) const;
ICoreMatrix elementWiseBetaInc(const ICoreMatrix& a, const ICoreMatrix& b,
const bool& upper, std::string* whyNot = nullptr) const;
// n x n / m x n filled with one value -- what Inf(m, n) and NaN(m, n) are
// built from.
static ICoreMatrix constant(const size_t& size1, const size_t& size2, const double& value);
// nchoosek(n, k) for a SCALAR n: the binomial coefficient -- MATLAB's other
// reading of the same name.
static double binomialCoefficient(const double& n, const double& k,
std::string* whyNot = nullptr);
// nchoosek(v, k): every k-element combination of v's entries, one per row,
// in MATLAB's lexicographic order. An empty result is 0 x k.
static ICoreMatrix combinations(const ICoreMatrix& vector, const size_t& k,
std::string* whyNot = nullptr);
// Element-wise comparisons (broadcast a scalar operand); result is 0.0/1.0 per entry
ICoreMatrix greaterThan(const ICoreMatrix& other) const;
ICoreMatrix lessThan(const ICoreMatrix& other) const;
ICoreMatrix greaterEqual(const ICoreMatrix& other) const;
ICoreMatrix lessEqual(const ICoreMatrix& other) const;
ICoreMatrix equalTo(const ICoreMatrix& other, double tol = 1e-12) const;
ICoreMatrix notEqualTo(const ICoreMatrix& other, double tol = 1e-12) const;
// Matrix power: integer exponent via exponentiation by squaring (negative uses inverse())
ICoreMatrix mpower(const long long& n) const;
// Reductions over all entries
double sum() const;
double mean() const;
double maxElement() const;
double minElement() const;
// Axis-wise reductions, 2D only. *Rows collapses each row to one value
// (result is a column vector); *Cols collapses each column (row vector).
ICoreMatrix sumRows() const;
ICoreMatrix sumCols() const;
ICoreMatrix meanRows() const;
ICoreMatrix meanCols() const;
ICoreMatrix maxRows() const;
ICoreMatrix maxCols() const;
ICoreMatrix minRows() const;
ICoreMatrix minCols() const;
// Cumulative sum/product. Vectors accumulate along their length; general
// matrices accumulate down each column (Matlab convention).
ICoreMatrix cumsum() const;
ICoreMatrix cumprod() const;
// Shape manipulation, 2D only. reshape reads and writes COLUMN-major,
// MATLAB's element order (C0.8) -- it was row-major until 2026-09-02.
ICoreMatrix reshape(const size_t& newDim1, const size_t& newDim2) const;
ICoreMatrix flipud() const;
ICoreMatrix fliplr() const;
ICoreMatrix repmat(const size_t& rowTimes, const size_t& colTimes) const;
// ---- Core MATLAB's array structure and shape (the core-MATLAB board, C4) --
//
// MATLAB semantics at the language level, whatever this class stores
// internally (the board's C0.8): 1-based counting is the evaluator's business, and the
// one ordering decision that reaches here is `nonzeros`, which is
// COLUMN-major because MATLAB's is.
//
// An operation whose MATLAB answer is an EMPTY matrix refuses instead,
// through whyNot, naming C3.9: the console has no empty-matrix semantics
// yet (`[] + 1` is `[]` in MATLAB and nothing here), so returning a 0 x n
// would hand the rest of the console a value it cannot compute with. Same
// contract the real* family uses -- a plausible answer would be the one
// answer that cannot be right.
// k-th diagonal, MATLAB's diag(v, k) / diag(A, k): a VECTOR builds an
// (n + |k|) square carrying v on that diagonal, a MATRIX returns that
// diagonal as a column. k = 0 is diag() above.
ICoreMatrix diag(const int& k, std::string* whyNot = nullptr) const;
// Block diagonal (MATLAB blkdiag): the two operands down the diagonal,
// zeros elsewhere. Neither need be square and their sizes need not match.
static ICoreMatrix blkdiag(const ICoreMatrix& first, const ICoreMatrix& second);
// Size, in MATLAB's reading: numel() counts every entry, length() is the
// longest dimension (0 when the matrix is empty), ndims() is 2 for every
// 2-D matrix -- MATLAB never reports fewer than 2, so a scalar is 2, not 0.
size_t numel() const;
size_t length() const;
size_t ndims() const;
// Upper/lower triangle at and beyond the k-th diagonal (MATLAB triu/tril).
ICoreMatrix triu(const int& k) const;
ICoreMatrix tril(const int& k) const;
// Counter-clockwise quarter turns (MATLAB rot90). `times` is taken mod 4
// and may be negative.
ICoreMatrix rot90(const int& times) const;
// Circular shift along `dim` (1 = down the columns, 2 = along the rows);
// `by` may be negative. MATLAB's no-dim form shifts along the first
// non-singleton dimension -- firstNonSingletonDim() below is that rule.
ICoreMatrix circshift(const long long& by, const int& dim) const;
// Reverse the order along `dim` (MATLAB flip); flipud and fliplr are its
// two fixed-dim spellings.
ICoreMatrix flip(const int& dim) const;
// n-th order difference along `dim` (MATLAB diff): each pass drops one
// entry from that dimension. An order that would consume the dimension
// refuses (C3.9) rather than returning an empty.
ICoreMatrix diff(const size_t& order, const int& dim,
std::string* whyNot = nullptr) const;
// The dimension MATLAB reduces along when none is named: the first whose
// extent is not 1, or 1 when there is none (a scalar).
int firstNonSingletonDim() const;
// Non-zero entries: nnz() counts them; nonzeros() returns them as a COLUMN
// in COLUMN-major order, MATLAB's order (C0.8). All-zero refuses (C3.9).
size_t nnz() const;
ICoreMatrix nonzeros(std::string* whyNot = nullptr) const;
// ---- Ordering, searching and sets (the core-MATLAB board, C4.24-C4.28,
// C9.5, C9.6) ------------------------------------------------------
//
// Three of MATLAB's rules here are the opposite of the obvious one, all
// measured against R2026a rather than reasoned about, and every one falls
// out of using plain `==` and an is-NaN test rather than a tolerance:
//
// * A NaN is never equal to itself, so `unique([NaN 1 NaN])` keeps BOTH
// NaNs (`[1 NaN NaN]`) and `ismember(NaN, NaN)` is FALSE.
// * A NaN sorts LAST ascending and FIRST descending -- it is not simply
// "largest", it moves with the direction.
// * Negative zero is zero: `unique([-0 0])` is one element, and
// `find([0 -0 1])` is 3.
//
// Ordering is column-major wherever an index or a flattening is involved
// (the board's C0.8), and an empty result REFUSES naming C3.9 rather than
// returning a 0 x n, as the C4.34 family does.
//
// Orientation follows MATLAB exactly and is not the same rule everywhere:
// `unique`/`find` give a ROW only for a row-vector input, while the set
// operations give a row only when BOTH inputs are rows.
// sort along `dim` (1 = down the columns, 2 = along the rows).
ICoreMatrix sorted(const int& dim, const bool& descending) const;
// sortrows: order whole ROWS lexicographically by the given 1-based
// columns, a negative column meaning descending on it. An empty `columns`
// means every column, left to right. Ties keep their original order.
ICoreMatrix sortedRows(const std::vector<long long>& columns,
std::string* whyNot = nullptr) const;
// unique: sorted-and-deduplicated (or first-appearance order when
// `stable`). `byRows` deduplicates whole rows instead of elements.
ICoreMatrix uniqueValues(const bool& stable, const bool& byRows,
std::string* whyNot = nullptr) const;
// find: the 1-based COLUMN-MAJOR indices of the non-zero entries, at most
// `limit` of them (0 = all), taken from the end when `fromLast`.
ICoreMatrix findIndices(const size_t& limit, const bool& fromLast,
std::string* whyNot = nullptr) const;
// sub2ind over a rows x cols shape: 1-based, column-major, element-wise
// over the two subscript matrices.
static ICoreMatrix subToInd(const size_t& rows, const size_t& cols,
const ICoreMatrix& rowSubs, const ICoreMatrix& colSubs,
std::string* whyNot = nullptr);
// ismember: 1.0 where the entry appears in `set`, shaped like *this.
// `byRows` asks the question of whole rows instead, one answer per row.
ICoreMatrix memberOf(const ICoreMatrix& set) const;
ICoreMatrix memberOfRows(const ICoreMatrix& set, std::string* whyNot = nullptr) const;
// The four set operations, over elements or over whole rows.
enum class SetOp { Union, Intersect, Difference, SymmetricDifference };
static ICoreMatrix setOperation(const ICoreMatrix& first, const ICoreMatrix& second,
const SetOp& op, const bool& stable,
std::string* whyNot = nullptr);
static ICoreMatrix setOperationRows(const ICoreMatrix& first, const ICoreMatrix& second,
const SetOp& op, std::string* whyNot = nullptr);
// ---- Core MATLAB's construction and scalar linear algebra (C4, C7) -----
//
// Same contract as the shape family above: MATLAB semantics at the
// language level, and an operation whose MATLAB answer is EMPTY refuses
// through whyNot rather than handing back a 0 x n the rest of the console
// cannot compute with (C3.9).
// Rectangular identity (MATLAB eye(m, n)): ones on the main diagonal.
static ICoreMatrix eye(const size_t& size1, const size_t& size2);
// MATLAB's magic(n) square, algorithm for algorithm: the Siamese method
// for odd n, the complement pattern when 4 divides n, and the LUX row
// exchange for the singly even n -- which is why magic(2) is MATLAB's
// degenerate [1 3; 4 2] rather than a refusal. n = 0 refuses (C3.9).
static ICoreMatrix magic(const size_t& n, std::string* whyNot = nullptr);
// MATLAB's colon(a, s, b), a ROW vector. The element COUNT is MATLAB's:
// an all-integer triple counts by truncation, anything else rounds and
// then steps back when a + n*s overshoots b by more than the tolerance.
// The last element is set to b EXACTLY when it lands inside that
// tolerance, which is why 0:0.1:1 ends at 1 and not at 0.9999999999999999.
// An empty range (s = 0, or b out of reach) refuses (C3.9).
static ICoreMatrix colon(const double& from, const double& step, const double& to,
std::string* whyNot = nullptr);
// Products, the mirror of sum()/sumRows()/sumCols(): over every entry,
// along each row (a column vector), down each column (a row vector).
double prod() const;
ICoreMatrix prodRows() const;
ICoreMatrix prodCols() const;
// MATLAB's norm(A, p). A VECTOR takes any p > 0 and both infinities
// (-Inf is the smallest magnitude); a MATRIX takes 1, 2 and Inf, since
// 'fro' is norm_Frobenius() above. Any other p refuses, as MATLAB errors.
double norm_p(const double& p, std::string* whyNot = nullptr) const;
// MATLAB's vecnorm: the p-norm of each column (dim = 1, a row vector) or
// of each row (dim = 2, a column vector).
ICoreMatrix vecnorm(const double& p, const int& dim, std::string* whyNot = nullptr) const;
// MATLAB's rcond, which is LAPACK's one-norm condition ESTIMATE and not
// the quantity it estimates -- so this is the estimator (DLACN2 driven the
// way DGECON drives it), not 1/(norm(A,1)*norm(inv(A),1)). The two differ
// by 15% on the parity suite's own 4x4 fixture, measured against R2026a
// 2026-09-02 (C7.7). A singular A answers 0, as MATLAB's does.
double rcond() const;
// MATLAB's cond(A, p) = norm(A, p) * norm(inv(A), p). p = 2 comes from the
// SVD (conditionNumber() above); p = -1 is this file's spelling of
// MATLAB's 'fro', the one p that is not a number.
double conditionNumber(const double& p, std::string* whyNot = nullptr) const;
// MATLAB's DEFAULT rank tolerance, max(m, n) * eps(largest singular
// value). rank(tol) above already counts the singular values above tol --
// MATLAB's rule; only the default differed, and this is that default.
double rankDefaultTolerance() const;
// MATLAB's inv(): a singular operand answers an all-Inf matrix and says so
// through wasSingular, where inverse() above returns the default matrix
// and logs an error. Both exist because this tree's own callers (the
// state-space discretizations, the descriptor block) test that empty.
ICoreMatrix inverseOrInf(bool* wasSingular = nullptr) const;
// ---- Cumulative reductions, integration and binning (C4.32, C4.37,
// ---- C4.38, C7.15, C8.16, C8.17) --------------------------------------
//
// Same contract as the two families above: MATLAB semantics at the
// language level, and an operation whose MATLAB answer is EMPTY refuses
// through whyNot rather than handing back a 0 x n (C3.9).
// Cumulative reductions along a dimension (MATLAB cumsum/cumprod/cummax/
// cummin, C4.32). `dim` is 1 (down the columns) or 2 (along the rows); the
// no-argument cumsum()/cumprod() above are the dim = firstNonSingletonDim()
// spellings, which is MATLAB's own default. Unlike the collapsing
// reductions these keep the shape of the input.
ICoreMatrix cumsum(const int& dim) const;
ICoreMatrix cumprod(const int& dim) const;
ICoreMatrix cummax(const int& dim) const;
ICoreMatrix cummin(const int& dim) const;
// Trapezoidal integration (MATLAB trapz/cumtrapz, C8.16/C8.17). trapz
// collapses `dim` to a single entry; cumtrapz keeps the shape and its
// first entry along `dim` is 0. The x-taking forms integrate over a
// non-uniform grid and need one x per sample along `dim` -- MATLAB does
// not require x to be sorted or increasing, and neither do these; the
// others use unit spacing.
ICoreMatrix trapz(const int& dim) const;
ICoreMatrix trapz(const ICoreMatrix& x, const int& dim,
std::string* whyNot = nullptr) const;
ICoreMatrix cumtrapz(const int& dim) const;
ICoreMatrix cumtrapz(const ICoreMatrix& x, const int& dim,
std::string* whyNot = nullptr) const;
// MATLAB's numeric gradient (C8.18): a CENTRAL difference at every
// interior sample and a ONE-SIDED difference at each end, so the answer
// has the shape of the input rather than one entry fewer (that is diff()).
// `dim` is 1 down the rows -- MATLAB's FY -- and 2 across the columns,
// MATLAB's FX and the one a single-output gradient(A) answers.
// `coords` is NULL for unit spacing, a SCALAR for uniform spacing, or one
// coordinate per sample along `dim`, in which case the interior difference
// is divided by the distance its two neighbours actually span. A lane of
// one sample is 0, as MATLAB's gradient(5) is.
ICoreMatrix gradient(const int& dim, const ICoreMatrix* coords,
std::string* whyNot = nullptr) const;
// Orthonormal basis for the range of A (MATLAB orth, C7.15): the first
// rank(A) left singular vectors, as columns. `tol` below zero asks for
// MATLAB's automatic threshold, max(m, n) * eps * sigma_max -- the same
// rule rankDefaultTolerance() above states. The basis is NOT unique (each
// column's sign, and the order inside a repeated singular value, are
// conventions), so a test pins the projector Q*Q', which is. A zero matrix
// has an empty range and refuses (C3.9).
ICoreMatrix orth(const double& tol = -1.0, std::string* whyNot = nullptr) const;
// accumarray(subs, val), MATLAB's basic two-argument form (C4.37): subs is
// an n x 1 column of 1-based row subscripts or an n x 2 of [row, col], and
// val is one value per subscript row or a scalar broadcast over them all.
// The result is max(subs) tall (and wide, for the two-column form), each
// entry the SUM of the values that landed on it and zero where none did.
static ICoreMatrix accumarray(const ICoreMatrix& subs, const ICoreMatrix& values,
std::string* whyNot = nullptr);
// histcounts (C4.38), in two halves so the binning rule can be tested on
// its own. histogramEdges reproduces MATLAB's OWN bin chooser -- the
// integer and Scott rules of toolbox/matlab/datafun/histcounts.m over
// .../+matlab/+internal/+math/binpicker.m -- so histcounts(x) and
// histcounts(x, nbins) are value-comparable against MATLAB rather than
// merely plausible; `bins` 0 asks for the automatic rule, anything else
// for that many bins. histcounts then counts x into `edges`: bin i is
// [e(i), e(i+1)), the LAST bin is closed, values outside every bin and NaN
// are not counted, and the result is a ROW one shorter than `edges`.
static ICoreMatrix histogramEdges(const ICoreMatrix& x, const size_t& bins,
std::string* whyNot = nullptr);
static ICoreMatrix histcounts(const ICoreMatrix& x, const ICoreMatrix& edges,
std::string* whyNot = nullptr);
// ---- Concatenation, dot products, correlation, spectral shift and the
// ---- two least-squares variants (C4.12, C5.38, C6.6, C7.19, C10.4) ----
// MATLAB's cat(dim, A, B): dim 1 stacks, dim 2 joins side by side -- which
// are concatVertical and concatHorizontal above, under MATLAB's spelling.
ICoreMatrix cat(const int& dim, const ICoreMatrix& other,
std::string* whyNot = nullptr) const;
// MATLAB's dot: the scalar product of two VECTORS, and column by column
// (dim 1) or row by row (dim 2) for two matrices. The dot() above is the
// console's older "over every element" reading and is NOT this one -- the
// two agree on vectors and disagree on everything else (the board's C0.1).
ICoreMatrix dotProduct(const ICoreMatrix& other, const int& dim,
std::string* whyNot = nullptr) const;
// MATLAB's corrcoef(x, y): the 2 x 2 correlation of two vectors read as
// two variables. correlation() above is corrcoef(A) for a data matrix,
// and this is the pair form, which is corrcoef([x(:) y(:)]).
static ICoreMatrix correlationOfPair(const ICoreMatrix& x, const ICoreMatrix& y,
std::string* whyNot = nullptr);
// MATLAB's fftshift/ifftshift: swap the halves of `dim` (0 = every
// dimension) so a zero-frequency term moves to the centre and back. Both
// are a circshift of floor(n/2) -- forward one way, backward the other --
// and that is the whole of it. They differ only for an ODD length, where
// the split is off centre and shifting the same way twice does not return.
ICoreMatrix fftshift(const int& dim, const bool& inverse) const;
// MATLAB's lsqminnorm(A, b): of all the least-squares solutions, the one
// of smallest norm. Equal to A\b whenever A has full column rank, and the
// pseudo-inverse solution when it does not -- which is exactly the case
// the backslash operator refuses to guess at (C1.3).
ICoreMatrix lsqminnorm(const ICoreMatrix& b, std::string* whyNot = nullptr) const;
// MATLAB's lscov(A, b) and lscov(A, b, w): ordinary and weighted least
// squares, the weighted form scaling every row of A and b by sqrt(w)
// first -- MATLAB's own reduction, from optimfun-adjacent lscov.m. Pass a
// 1 x 1 of 1.0 for the unweighted form. A rank-deficient design refuses:
// MATLAB warns and returns a BASIC solution, which is a different vector
// from the minimum-norm one and nothing in the call says which was meant.
ICoreMatrix lscov(const ICoreMatrix& b, const ICoreMatrix& weights,
std::string* whyNot = nullptr) const;
// ---- The one random generator, and MATLAB's rng (C4.4-C4.7) -----------
//
// rand, randn, randi and randperm all draw from ONE Mersenne Twister,
// because rng(seed) has to mean something for all four at once. It starts
// seeded 0, which is MATLAB's own startup state (rng('default')), so a
// fresh console replays a fresh MATLAB's SEQUENCE OF DRAWS -- not its
// values: MATLAB's stream and mt19937's never agree sample for sample,
// which is why every parity row in this family is Skip (C0.9). What the
// seed buys, and what the regress cases assert, is that the same seed
// replays the same stream.
static unsigned int randomSeed(); // the seed now in force
static unsigned int seedRandom(const unsigned int& seed); // set it; returns the PREVIOUS seed
static unsigned int shuffleRandom(); // seed unpredictably; returns the PREVIOUS seed
// MATLAB's randi: whole numbers drawn uniformly from the CLOSED range
// [low, high]. Both ends must be whole numbers and low must not exceed
// high, as MATLAB requires.
static ICoreMatrix randi(const double& low, const double& high,
const size_t& size1, const size_t& size2,
std::string* whyNot = nullptr);
// MATLAB's randperm(n) / randperm(n, k): k values drawn from 1..n WITHOUT
// replacement, as a ROW -- MATLAB's shape, and a new name takes it from
// birth. k = n is the full permutation. k > n refuses, as MATLAB errors.
static ICoreMatrix randperm(const size_t& n, const size_t& k,
std::string* whyNot = nullptr);
// MATLAB's lsqnonneg(C, d): the least-squares solution of C*x = d subject
// to x >= 0, by the Lawson-Hanson active-set algorithm transcribed from
// optimfun/lsqnonneg.m -- its convergence tolerance included,
// 10 * eps * norm(C, 1) * max(size(C)), because that tolerance is what
// decides which variables come back exactly zero (C7.20).
ICoreMatrix lsqnonneg(const ICoreMatrix& d, std::string* whyNot = nullptr) const;
// Convolution (MATLAB conv and conv2). `shape` is "full" (the whole
// result), "same" (the central part, the size of the FIRST operand) or
// "valid" (only the positions where the operands fully overlap). A
// "valid" result that would be empty refuses (C3.9), as does a 1-D
// convolution of anything but two vectors, which MATLAB refuses too.
//
// MATLAB's orientation rule for the 1-D form: the answer is a column only
// when BOTH operands are column-shaped; a row anywhere makes it a row.
static ICoreMatrix convolve(const ICoreMatrix& first, const ICoreMatrix& second,
const std::string& shape, std::string* whyNot = nullptr);
static ICoreMatrix convolve2D(const ICoreMatrix& first, const ICoreMatrix& second,
const std::string& shape, std::string* whyNot = nullptr);
// MATLAB's filter (the core-MATLAB row C10.6): the rational difference
// equation
//
// a(1)y(n) = b(1)x(n) + ... + b(N)x(n-N+1) - a(2)y(n-1) - ... - a(N)y(n-N+1)
//
// run in the DIRECT FORM II TRANSPOSED structure MATLAB itself uses. The
// structure is not an implementation detail here: it is what defines the
// state vector, so zi and zf mean what MATLAB's mean only under this form.
// b and a are vectors, zero-padded to the longer of the two; a(1) must be
// non-zero and both are divided through by it, as MATLAB does.
//
// `dim` is 1 (down each column) or 2 (along each row); 0 asks for MATLAB's
// default, the first non-singleton dimension -- 2 for a row vector, 1 for
// everything else.
//
// `initial` is MATLAB's zi (nullptr for none) and `finalState` its zf
// (nullptr when the caller does not want it). Both carry the state along
// the LEADING dimension whatever `dim` is: with N-1 states and S signals
// they are (N-1) BY S, even when dim is 2 and the signals run across the
// rows -- so filtering a 3x4 array along dim 2 answers a 2x3 zf, not 3x2.
// Measured against MATLAB R2026a on 2026-09-02; MATLAB's own error message
// for a mis-shaped zi says it in as many words ("an array with the leading
// dimension of size max(length(a),length(b))-1 and with remaining
// dimensions matching those of x").
//
// A single signal takes any vector of N-1 initial states, whichever way it
// is turned; several signals need the (N-1) BY S array exactly.
static ICoreMatrix filter(const ICoreMatrix& numerator, const ICoreMatrix& denominator,
const ICoreMatrix& x, const ICoreMatrix* initial,
const int& dim, ICoreMatrix* finalState,
std::string* whyNot = nullptr);
// ---- The discrete field operators (the core board's C8.27) ------------
//
// All three are core MATLAB -- `which del2` is datafun/, `which curl` and
// `which divergence` are graphics/specgraph/ -- so the board's "🚫 N-D
// fields" was a statement about their THREE-dimensional forms only. The
// 2-D forms are these, and the N-D ones stay out under Z6.
// MATLAB's del2: the discrete Laplacian, DIVIDED BY 2*ndims. It is
// transcribed from `toolbox/matlab/datafun/del2.m` rather than derived,
// because two of its rules are not what a derivation would produce:
//
// * the divisor is `ndims(f)` AFTER MATLAB has made a vector into a
// column, and MATLAB's ndims is never below 2 -- so a vector is
// divided by 2 as a matrix is, and del2([1 4 9 16]) is 0.5 and not 1
// even though only one dimension was differenced.
// * the end samples are a LINEAR EXTRAPOLATION of the interior second
// differences when there are more than three samples along the
// dimension, a copy of the single interior value when there are
// exactly three, and 0 when there are two or fewer.
//
// `hx`/`hy` are NULL for unit spacing, a scalar step, or one coordinate
// per sample; del2(f, h) applies the one step to BOTH dimensions, which is
// what MATLAB's scalar form does. A row vector answers a row.
static ICoreMatrix laplacian(const ICoreMatrix& f, const ICoreMatrix* hx,
const ICoreMatrix* hy, std::string* whyNot = nullptr);
// MATLAB's 2-D divergence: d(u)/dx + d(v)/dy, each term the numeric
// gradient above. `x`/`y` are the meshgrid pair (NULL for unit spacing);
// only x's first ROW and y's first COLUMN are read, which is what
// MATLAB's own divergence.m does with them.
static ICoreMatrix divergence2D(const ICoreMatrix& u, const ICoreMatrix& v,
const ICoreMatrix* x, const ICoreMatrix* y,
std::string* whyNot = nullptr);
// MATLAB's 2-D curl: d(v)/dx - d(u)/dy.
//
// **This function answers a DIFFERENT number depending on how many
// outputs are asked for**, and it is MATLAB's rule rather than a
// convenience: `curl(X, Y, U, V)` with ONE output is the ANGULAR VELOCITY
// -- half the curl -- while `[curlz, cav] = curl(X, Y, U, V)` answers the
// curl first and the angular velocity second. Read off curl.m and
// measured against R2026a on 2026-09-02: for U = X.*Y, V = X.^2 on
// meshgrid(1:4, 1:4) the one-output form is 1, 1, 1.5, 1.5 across each
// row where the curl itself is 2, 2, 3, 3. `angularVelocityOnly` says
// which the caller is asking for.
static ICoreMatrix curl2D(const ICoreMatrix& u, const ICoreMatrix& v,
const ICoreMatrix* x, const ICoreMatrix* y,
const bool& angularVelocityOnly,
std::string* whyNot = nullptr);
// ---- The two signal functions that are core MATLAB (C10.10) -----------
//
// `which detrend` and `which xcorr` are BOTH toolbox/matlab/datafun/ in
// R2026a -- base MATLAB, not the Signal Processing Toolbox -- which is
// the measurement the board's row asked for before marking itself 🚫.
// The other seven names on that row are signal/signal/ and stay refused.
// MATLAB's detrend: subtract the least-squares polynomial of degree
// `degree` from each COLUMN (from the vector, whichever way it is
// turned). Degree 0 is MATLAB's 'constant' and is the column mean
// subtracted exactly; degree 1 is the default and MATLAB's 'linear'.
//
// The fit is over the SCALED sample points detrend.m uses -- s./s(end),
// so the last one is 1 -- rather than over 1..N. It is the same fit in
// exact arithmetic and a better-conditioned one in this one, and matching
// MATLAB's conditioning is the difference between agreeing to the last
// bit and agreeing to the tenth digit on a long ramp.
static ICoreMatrix detrend(const ICoreMatrix& x, const size_t& degree,
std::string* whyNot = nullptr);
// MATLAB's xcorr, for VECTORS: the raw cross-correlation
//
// c(lag) = sum_n x(n + lag) * y(n)
//
// at every lag from -maxlag to +maxlag, so the answer has 2*maxlag+1
// entries and c's middle entry is lag 0. `y` empty asks for the
// AUTOcorrelation. `maxlag` below zero asks for MATLAB's default,
// max(numel(x), numel(y)) - 1; a larger one zero-pads, as MATLAB's does.
//
// `scale` is "none" (the default), "biased" (divided by N), "unbiased"
// (divided by N - |lag|, with a non-positive divisor clamped to 1, which
// is MATLAB's own guard) or "coeff"/"normalized" -- the autocorrelation
// divided by its own zero-lag value, or the cross-correlation divided by
// sqrt(sum(x.^2) * sum(y.^2)). Every scale but "none" needs the two
// vectors the same length, which is MATLAB's NoScale error.
//
// The orientation of the ANSWER follows x: a row answers a row, a column a
// column. The LAGS are a row either way, which is MATLAB's shape and the
// one asymmetry here -- a lag axis is not a second signal.
//
// This is the direct O(N * lags) sum where MATLAB's is an FFT and an
// inverse FFT, so the two agree to rounding rather than exactly -- and it
// is MATLAB that carries the error: xcorr([1 2 3 4]) answers
// 3.9999999999999982 at its end lags where the sum of products is 4.
// `y` is NULL for the autocorrelation. It is a pointer rather than an
// empty matrix because a default-constructed ICoreMatrix is a 1 x 1
// holding ZERO, not an empty one -- the same trap that made chol() answer
// a silent 0 before C7.12.
static ICoreMatrix xcorr(const ICoreMatrix& x, const ICoreMatrix* y,
const int& maxlag, const std::string& scale,
ICoreMatrix* lags, std::string* whyNot = nullptr);
// Statistics: columns are variables, rows are observations
ICoreMatrix covariance() const;
ICoreMatrix correlation() const;
// Vector-only helpers
ICoreMatrix normalize() const; // unit vector (v / ||v||)
double angleTo(const ICoreMatrix& other) const; // angle between two vectors, radians
// Sylvester-family linear matrix equations, solved via vec()/Kronecker
// reduction. A unique solution needs the operator nonsingular — lyap:
// eig(A)+eig(A) never 0; dlyap: eig(A)*eig(A) never 1; sylvester:
// eig(A)+eig(B) never 0. A singular operator makes the QR solve return a
// NON-solution with no diagnostic, so each solver verifies its residual
// afterwards; on failure it reports (into whyNot when given, the run
// diagnosis otherwise) and returns the default matrix. Found by the
// command-parity suite: MATLAB returns NaN on the same inputs.
ICoreMatrix lyap(const ICoreMatrix& Q, std::string* whyNot = nullptr) const; // A*X + X*A^T + Q = 0
ICoreMatrix dlyap(const ICoreMatrix& Q, std::string* whyNot = nullptr) const; // A*X*A^T - X + Q = 0
ICoreMatrix sylvester(const ICoreMatrix& B, const ICoreMatrix& C,
std::string* whyNot = nullptr) const; // A*X + X*B = C
// The GENERALIZED (descriptor) members of the same family, and the
// discrete Sylvester form -- MATLAB's lyap(A, Q, [], E), dlyap(A, Q, [], E)
// and dlyap(A, B, C) (the toolbox board's T1.28). Same vec()/
// Kronecker reduction, same rank-and-residual verification, so a singular
// operator is reported rather than answered:
// lyapGeneralized A*X*E' + E*X*A' + Q = 0 (E⊗A + A⊗E) vec(X) = -vec(Q)
// dlyapGeneralized A*X*A' - E*X*E' + Q = 0 (A⊗A - E⊗E) vec(X) = -vec(Q)
// dlyapSylvester A*X*B - X + C = 0 (B'⊗A - I) vec(X) = -vec(C)
//
// Note the SIGNS, which are MATLAB's and not this family's: the two
// Lyapunov forms take Q on the left-hand side and answer the X that makes
// the sum ZERO, while sylvester() above takes C on the RIGHT and answers
// A*X + X*B = C. That is why MATLAB's lyap(A, B, C) is this console's
// sylvester(A, B, -C) and not a fourth solver.
ICoreMatrix lyapGeneralized(const ICoreMatrix& Q, const ICoreMatrix& E,
std::string* whyNot = nullptr) const;
ICoreMatrix dlyapGeneralized(const ICoreMatrix& Q, const ICoreMatrix& E,
std::string* whyNot = nullptr) const;
ICoreMatrix dlyapSylvester(const ICoreMatrix& B, const ICoreMatrix& C,
std::string* whyNot = nullptr) const;
// Algebraic Riccati equations (LQR/Kalman-style design). Solved via the stable
// invariant subspace of the Hamiltonian (care) / symplectic pencil (dare).
// Assumes (A,B) stabilizable and (A,Q) detectable, with no eigenvalues exactly
// on the stability boundary (imaginary axis / unit circle).
ICoreMatrix care(const ICoreMatrix& B, const ICoreMatrix& Q, const ICoreMatrix& R) const; // A^T X + X A - X B R^-1 B^T X + Q = 0
ICoreMatrix dare(const ICoreMatrix& B, const ICoreMatrix& Q, const ICoreMatrix& R) const; // A^T X A - X - A^T X B (R+B^T X B)^-1 B^T X A + Q = 0
Eigen::MatrixXd toEigen2DMatrixXd() const;
////////////////////////////////////////////////////
///
/// Matrix print-out operator
///
////////////////////////////////////////////////////
// Declare friend operator<< for output streaming
friend std::ostream& operator<<(std::ostream& os, const ICoreMatrix& mat);
// A matrix is copied constantly -- by value into and out of nearly every
// numeric call in the tree -- so the type stays copyable and nothrow-movable.
// The unique_ptr residue below deletes the implicit copy, so all four are
// written out in the .cpp.
//
// THE COST, recorded rather than hidden: an ICoreMatrix now costs TWO heap
// allocations to build instead of one (the Impl, then its data vector) and
// one more indirection per element access. It already allocated for `data`
// on every copy, so this is a constant factor on an allocating path, not a
// new allocation on a free one. See the SH2 notes in HEADER_SURFACE.md.
ICoreMatrix(const ICoreMatrix& other);
ICoreMatrix& operator=(const ICoreMatrix& other);
ICoreMatrix(ICoreMatrix&& other) noexcept;
ICoreMatrix& operator=(ICoreMatrix&& other) noexcept;
// Declared, defined in the .cpp: unique_ptr cannot destroy an incomplete Impl.
~ICoreMatrix();
private:
class Impl; // the two-line residue; state lives here
std::unique_ptr<Impl> impl;
};ICorePolynomial.h#
src/ICoreBlocks/ICoreMath/Foundation/ICorePolynomial.h
MATLAB's polyval and polyvalm (the core-MATLAB board, C8.2 and C8.10). polyval evaluates at EVERY entry of x and keeps its shape; polyvalm evaluates the polynomial IN the matrix -- x^2 is A*A, not A .* A, and the constant term is that multiple of the identity, so A must be square.
ICorePolynomial#
ICorePolynomial.h:9 · class · 25 declaration(s)
class ICorePolynomial {
public:
ICorePolynomial();
explicit ICorePolynomial(const double& coefficient);
explicit ICorePolynomial(const std::vector<double>& coefficients);
explicit ICorePolynomial(const ICoreMatrix& coefficients);
double evaluate(const double& x) const;
// MATLAB's polyval and polyvalm (the core-MATLAB board, C8.2 and C8.10).
// polyval evaluates at EVERY entry of x and keeps its shape; polyvalm
// evaluates the polynomial IN the matrix -- x^2 is A*A, not A .* A, and
// the constant term is that multiple of the identity, so A must be square.
ICoreMatrix evaluateElementWise(const ICoreMatrix& x) const;
ICoreMatrix evaluateMatrix(const ICoreMatrix& matrix, std::string* whyNot = nullptr) const;
ICoreComplexVariable evaluateComplex(const double &real, const double &imag) const;
ICorePolynomial derivative() const;
ICorePolynomial integral(const double &constant) const;
void normalize();
void trimTrailingZeros(double tol = 1e-14);
ICorePolynomial operator+(const ICorePolynomial &other) const;
ICorePolynomial operator-(const ICorePolynomial &other) const;
ICorePolynomial operator*(const ICorePolynomial &other) const;
std::pair<ICorePolynomial, ICorePolynomial> divide(const ICorePolynomial &divisor) const;
// MATLAB's polyfit (the core-MATLAB board, C8.5): the degree-n polynomial
// whose coefficients least-squares fit the samples, answered as a ROW of
// n+1 coefficients highest power first -- MATLAB has no polynomial object
// to hand back, and a column x and y still answer a row.
//
// The fit is the backslash solution of the Vandermonde system, which is
// what MATLAB's own polyfit.m computes (a column-pivoting QR).
//
// `degree` must be strictly below the number of samples. MATLAB warns and
// answers a rank-deficient particular solution instead; that solution is
// not the fit anyone meant, so this refuses with MATLAB's own condition.
static ICoreMatrix fit(const ICoreMatrix& x, const ICoreMatrix& y, const size_t& degree,
std::string* whyNot = nullptr);
// MATLAB's residue (the core-MATLAB row C8.11): the partial fraction
// expansion of b(s)/a(s),
//
// b(s)/a(s) = r(1)/(s - p(1)) + ... + r(n)/(s - p(n)) + k(s)
//
// with a repeated pole taking consecutive terms of RISING power --
// r(j)/(s - p)^1, r(j+1)/(s - p)^2, ... -- which is MATLAB's own layout
// and is why r and p are the same length even when p repeats.
//
// The residues are solved as a LINEAR SYSTEM rather than through the
// derivative formula: multiplying the expansion by a(s) turns every term
// into the polynomial a(s)/(s - p(j))^e(j), so the coefficients of b give
// n equations in the n unknown residues. That is what MATLAB's own
// residue.m does, and it is the reason a repeated pole and a complex one
// need no special case here -- the derivative formula needs both.
//
// `poles` is roots(a) with equal poles adjacent; `residues` follows it
// term for term; `direct` is the polynomial part, EMPTY when b is of lower
// degree than a (MATLAB answers a 0x0 there, and so does this since C3.9).
// All three come back as columns, MATLAB's shape.
//
// Poles are complex in general, so r and p are complex too and narrow to
// real at the console boundary when nothing imaginary survives (C0.2).
static bool residue(const ICorePolynomial& numerator, const ICorePolynomial& denominator,
ICoreComplexMatrix& residues, ICoreComplexMatrix& poles,
ICoreMatrix& direct, std::string* whyNot = nullptr);
static ICorePolynomial fromRoots(const std::vector<ICoreComplexVariable> &roots);
std::vector<ICoreComplexVariable> roots() const;
ICoreMatrix companionMatrix() const;
bool isStableContinuous() const;
bool isStableDiscrete() const;
ICoreMatrix getCoefficients() const;
int getOrder() const;
void print(const std::string& varName = "x") const;
~ICorePolynomial();
};ICoreRecord.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreRecord.h
ICoreRecord#
ICoreRecord.h:45 · class · pImpl · 23 declaration(s)
MATLAB's scalar struct, as much of one as this console admits: an ORDERED, NAMED list of numeric fields, read by .name and nothing else (the toolbox board's T0.6 (A) and X3).
class ICoreRecord {
public:
// NOT `= default`: the residue below is a unique_ptr, and a defaulted
// default constructor leaves it null -- which compiles, links, and
// dereferences null on first use.
ICoreRecord();
// Records are stored in the variables space and handed around by value, so
// all four copy/move operations are written out in the .cpp.
ICoreRecord(const ICoreRecord& other);
ICoreRecord& operator=(const ICoreRecord& other);
ICoreRecord(ICoreRecord&& other) noexcept;
ICoreRecord& operator=(ICoreRecord&& other) noexcept;
~ICoreRecord();
void clear();
// Appends in the order MATLAB answers the fields. A name already present is
// REPLACED in place rather than appended twice, so a caller that builds a
// record in a loop cannot produce two fields with one name -- MATLAB's
// struct cannot hold that and neither can the storage string.
void append(const std::string& name, const ICoreMatrix& value);
void appendScalar(const std::string& name, const double& value);
// Records a field MATLAB answers that this console does not hold, so that
// reading it refuses BY NAME. `why` is the whole diagnostic the reader
// prints; it should name the rule and the row (see the class note).
void withhold(const std::string& name, const std::string& why);
[[nodiscard]] bool isEmpty() const;
[[nodiscard]] size_t fieldCount() const;
[[nodiscard]] const std::string& fieldName(const size_t& index) const;
[[nodiscard]] const ICoreMatrix& fieldValue(const size_t& index) const;
[[nodiscard]] std::vector<std::string> fieldNames() const;
[[nodiscard]] bool hasField(const std::string& name) const;
[[nodiscard]] bool getField(const std::string& name, ICoreMatrix& out) const;
// True when `name` is a field MATLAB answers here but this console withheld;
// `why` comes back with the diagnostic to print. Kept separate from
// hasField() so a reader cannot confuse "absent" with "refused".
[[nodiscard]] bool isWithheld(const std::string& name, std::string& why) const;
[[nodiscard]] std::vector<std::string> withheldNames() const;
// MATLAB's scalar-struct display: names right-aligned to the longest, one
// field per line, a non-scalar shown as its SIZE and not its values --
// `lower: [2x1 double]` -- exactly as MATLAB does. The numbers are
// MATLAB's `format short`, which this console already reproduces wherever
// a display is MATLAB's (a tf prints 0.3333, not 0.333333).
//
// ⚠ This display shows the numeric fields ONLY, so for a solver record
// carrying withheld character fields it is SHORTER than MATLAB's. That is
// deliberate and it is why the parity cases for these rows compare
// `out.iterations` and not `disp(out)` -- printing `algorithm: '...'` for
// a field the console refuses to hand back would be a display promising a
// value that no read can reach.
[[nodiscard]] std::string matlabDisplay() const;
// The lossless storage/bridge form: `struct('a', 1, 'b', [1, 2])`, with
// every number written at full precision. This is what the variables space
// stores and what crosses to MATLAB (X3); the import already admits a
// struct(...) call for the time-series shape, so the spelling is not new.
[[nodiscard]] std::string toStructCall() const;
static bool fromStructCall(const std::string& text, ICoreRecord& out, std::string& error);
// MATLAB's `format short` for ONE scalar, which is the display rule for a
// struct field and is NOT the four-significant-digit rule a tf's
// coefficients use. Public because the record is not the only caller that
// will want it, and two implementations of one rule is how they drift.
static std::string formatShort(const double& v);
private:
class Impl; // the two-line residue; state lives here
std::unique_ptr<Impl> impl;
};ICoreSpline.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreSpline.h
ICoreSpline#
ICoreSpline.h:12 · class · nested PiecewiseCubic · 0 declaration(s)
Interpolation through (x,y) sample points, backed by Eigen::Spline / Eigen::SplineFitting.
class ICoreSpline {
public:
// The B-spline fit behind the Studio Spline Fitting tool -- a spline of a
// CHOSEN degree through the samples, which is the one thing on this class
// MATLAB has no name for. It was called interp1() until C8.12 gave that
// name MATLAB's meaning; the degree argument it reads is what made the two
// impossible to keep under one name, since MATLAB's fourth argument is a
// method WORD.
//
// x, y: vectors of the same length (>= 2), x need not be evenly spaced but
// must be strictly increasing. xq: query points. degree: spline degree
// (default 3 = cubic); reduced automatically if there are too few points.
// Returns yq, one value per entry of xq. Queries outside [x(1), x(end)]
// are CLAMPED to the ends rather than extrapolated.
static ICoreMatrix bSplineFit(const ICoreMatrix& x, const ICoreMatrix& y,
const ICoreMatrix& xq, size_t degree = 3);
// MATLAB's interp1 (the core-MATLAB board, C8.12), whose default is
// LINEAR -- the divergence C0.1 decided in MATLAB's favour.
//
// The six methods are the ones with a single meaning. MATLAB's 'cubic' is
// deliberately absent: measured against R2026a 2026-09-02 it is a THIRD
// curve, neither Spline nor Pchip (on the uniform grid x = 1:4,
// v = [1 3 2 5] the three answer 2.375, 2.8125 and 2.4375 at xq = 1.5),
// and on a non-uniform grid it warns and silently becomes 'spline' while
// still refusing to extrapolate. Mapping it onto either of these would be
// a wrong number rather than a missing feature, so it is refused by name.
//
// `extrapolate` is MATLAB's 'extrap': false fills every query outside
// [min(x), max(x)] with `outsideValue` (MATLAB's extrapval, NaN when the
// caller passes none). MATLAB's DEFAULT differs by method -- Spline and
// Pchip extrapolate, the other four do not -- and that choice belongs to
// the caller, which is why it is a parameter here rather than a rule.
// Two measured asymmetries survive `extrapolate`: Previous has no previous
// sample BELOW the range and Next has no next one above, so those stay NaN
// however the flag is set.
//
// x and v are vectors of the same length (>= 2) with distinct entries;
// neither needs to be sorted (this sorts on x, as MATLAB does) and a
// decreasing x is legal. The answer has the SHAPE of xq.
enum class Method { Linear, Nearest, Next, Previous, Spline, Pchip };
static ICoreMatrix interp1(const ICoreMatrix& x, const ICoreMatrix& v,
const ICoreMatrix& xq, const Method& method,
const bool& extrapolate, const double& outsideValue,
std::string* whyNot = nullptr);
// MATLAB's spline() and pchip() (the core-MATLAB board, C8.14). Neither
// goes through Eigen: both are the piecewise cubic MATLAB itself builds,
// and they differ only in how the slope at each sample point is chosen.
//
// spline the NOT-A-KNOT cubic: the third derivative is continuous
// across the second and the second-to-last points, which is what
// makes two points a straight line and three a single parabola.
// pchip the shape-preserving cubic: the slope at a sample is a
// weighted harmonic mean of its neighbouring secants, and ZERO
// wherever they change sign, so the result never overshoots.
//
// Both EXTRAPOLATE outside [x(1), x(end)] by continuing the end cubic,
// as MATLAB does. x need not be sorted (both sort it, as MATLAB does) but
// must hold no repeats; x and y must be the same length.
static ICoreMatrix splineNotAKnot(const ICoreMatrix& x, const ICoreMatrix& y,
const ICoreMatrix& xq, std::string* whyNot = nullptr);
static ICoreMatrix pchip(const ICoreMatrix& x, const ICoreMatrix& y,
const ICoreMatrix& xq, std::string* whyNot = nullptr);
// The same two curves handed over as a FORM rather than evaluated: the
// sorted samples and the slope at each of them, which with the Hermite
// basis is the piecewise cubic itself (toolbox row T8.4). `spline` and
// `pchip` above are these two plus an evaluation, so there is one
// implementation of each curve and not two -- a caller that needs the
// COEFFICIENTS (to differentiate a fit exactly, or to integrate one)
// builds them from here rather than re-deriving the curve.
struct PiecewiseCubic {
std::vector<double> x; // sorted, no repeats
std::vector<double> y;
std::vector<double> slopes; // one per sample
};
static bool notAKnotForm(const ICoreMatrix& x, const ICoreMatrix& y,
PiecewiseCubic& out, std::string* whyNot = nullptr);
static bool pchipForm(const ICoreMatrix& x, const ICoreMatrix& y,
PiecewiseCubic& out, std::string* whyNot = nullptr);
// MATLAB's interp2 (the core-MATLAB board, C8.13) on a plaid grid.
//
// It is the SEPARABLE application of interp1 above, not a second
// algorithm: interpolate every grid column over y at the query rows, then
// interpolate those across x at the query columns. For each of the four
// methods that is exactly what MATLAB's own 2-D interpolant is -- a
// tensor product factors that way by construction -- so the two spellings
// cannot drift, the way conv2 and filter2 share one convolution.
//
// X and Y are either the meshgrid arrays (every row of X the same, every
// column of Y the same -- MATLAB refuses a grid that is not plaid, and so
// does this) or the two vectors they were made from. V is size(Y) BY
// size(X). Xq and Yq broadcast against each other the way MATLAB's do: the
// same shape is element-wise, a row against a column is the grid of pairs.
//
// The four methods are the ones with a single meaning in TWO dimensions,
// and the list is shorter than interp1's for reasons measured against
// R2026a on 2026-09-02 rather than assumed:
//
// previous, pchip MATLAB accepts the word, WARNS that it is 1-D only,
// and silently interpolates linearly instead. A word
// that changes method behind the caller's back is the
// one answer that cannot be right, so it is refused.
// cubic, v5cubic a THIRD curve again (C8.12): on the uniform grid
// x = y = 1:4 the three answer 1.265625 (cubic),
// 0.785156 (spline) and 1.334248 (makima) at
// (1.5, 1.5). On a non-uniform grid it also warns and
// becomes spline.
// makima a fourth curve, which interp1 does not have either.
//
// `extrapolate` and `outsideValue` mean what they mean on interp1, and
// MATLAB's default is per method there too: Spline extrapolates, the other
// three answer `outsideValue` (NaN unless the caller passes one) at any
// query outside the grid in EITHER coordinate.
enum class Method2D { Linear, Nearest, Next, Spline };
static ICoreMatrix interp2(const ICoreMatrix& gridX, const ICoreMatrix& gridY,
const ICoreMatrix& values,
const ICoreMatrix& queryX, const ICoreMatrix& queryY,
const Method2D& method, const bool& extrapolate,
const double& outsideValue, std::string* whyNot = nullptr);
};
};ICoreStatistics.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreStatistics.h
ICoreStatistics#
ICoreStatistics.h:24 · class · 1 declaration(s)
Core MATLAB's descriptive statistics (the core-MATLAB board, C6): median, var/std with MATLAB's weight argument, cov's three readings, rescale, and the mov* sliding-window family.
class ICoreStatistics {
public:
// The four reductions MATLAB spells with one name and an axis (C5.32-C5.35):
// sum, mean, max and min. `dim` is 1, 2 or 0 for the default (the first
// non-singleton dimension), `everyElement` is MATLAB's "all", and
// `omitNaN` is its "omitnan".
//
// **The NaN rules differ by function and they are MATLAB's, not a choice
// made here**: sum and mean PROPAGATE a NaN (one NaN makes the whole slice
// NaN) unless "omitnan" is asked for, while max and min IGNORE one always
// -- `max([1 NaN 3])` is 3, and "omitnan" only says out loud what max
// already does. A slice that is ALL NaN answers NaN whichever way.
//
// `indices`, when non-null, receives MATLAB's second output for max/min:
// the 1-based position of the chosen entry WITHIN its slice, or the linear
// (column-major) position over the whole array when `everyElement` is set.
// It is left untouched for sum and mean, which have no second output.
enum class Reduction { Sum, Mean, Max, Min };
static ICoreMatrix reduce(const Reduction& which, const ICoreMatrix& a, int dim,
bool everyElement, bool omitNaN,
ICoreMatrix* indices = nullptr, std::string* whyNot = nullptr);
// median(A), median(A, dim), median(..., "omitnan"). A NaN anywhere in a
// slice makes that slice's answer NaN unless omitNaN is set, and a slice
// that is ALL NaN answers NaN either way.
static ICoreMatrix median(const ICoreMatrix& a, int dim, bool omitNaN,
std::string* whyNot = nullptr);
// var(A), var(A, w), var(A, w, dim) -- and std, which is its square root
// exactly as MATLAB computes it. `weights` is null for MATLAB's w = 0
// (normalize by N-1), a 1x1 holding 0 or 1 for the flag, or a vector as
// long as the reduced dimension for the weighted variance (normalized by
// sum(w)). A one-element slice answers 0, not NaN, which is why var(5) is
// 0 in MATLAB and here.
static ICoreMatrix variance(const ICoreMatrix& a, const ICoreMatrix* weights,
int dim, std::string* whyNot = nullptr);
static ICoreMatrix standardDeviation(const ICoreMatrix& a, const ICoreMatrix* weights,
int dim, std::string* whyNot = nullptr);
// cov(A): a matrix is observations x variables and answers the variable
// covariance matrix; a VECTOR is one variable and answers its variance.
static ICoreMatrix covariance(const ICoreMatrix& a, std::string* whyNot = nullptr);
// cov(A, w) and cov(x, y) are the same call in MATLAB, told apart by the
// second argument: a scalar 0 or 1 is the normalization flag, anything
// else is a second variable, flattened column-major and paired with the
// first into the 2x2 covariance matrix.
static ICoreMatrix covariance(const ICoreMatrix& a, const ICoreMatrix& second,
std::string* whyNot = nullptr);
// mode(A), mode(A, dim), and the COUNT that [M, F] = mode(A) returns
// beside it (the core-MATLAB board, C6.2). Ties go to the SMALLEST value,
// as in MATLAB, and a NaN is not a candidate however often it repeats --
// mode([1, NaN, NaN, 2]) is 1. `counts` may be null when the caller wants
// only the value.
static ICoreMatrix mode(const ICoreMatrix& a, int dim, ICoreMatrix* counts,
std::string* whyNot = nullptr);
// bounds(A) -> [lo, hi] (C5.36): the column-wise minimum and maximum in
// one pass, or the whole array's when `everyElement` is set ("all"). NaN is
// omitted the way min/max omit it, and a slice that is all NaN answers NaN.
static bool boundsOf(const ICoreMatrix& a, int dim, bool everyElement,
ICoreMatrix& low, ICoreMatrix& high, std::string* whyNot = nullptr);
// rescale(A), rescale(A, l, u): the whole array is mapped onto [l, u] from
// its own min and max. Two MATLAB rules that are not obvious: a constant
// array answers l everywhere (not NaN), and an array holding an infinity
// answers NaN everywhere.
static ICoreMatrix rescale(const ICoreMatrix& a, const double& lo, const double& hi,
std::string* whyNot = nullptr);
// MATLAB's normalize (the core-MATLAB board, C6.8), which is NOT the
// console's older normalize() -- that one was a unit vector and now
// answers to `unit` (C0.1, decision A).
//
// ZScore (x - mean) / std, MATLAB's default and the N-1 std
// Center x - mean
// Scale x / std
// Range onto [lo, hi] from the slice's own min and max
// Norm x / norm(slice, p)
//
// `parameter` carries the method's own argument -- the p of Norm (default
// 2, and Inf is the largest magnitude) or the two-entry [lo hi] of Range
// (default [0 1]) -- and is null for the three that take none.
//
// **Every statistic here OMITS NaN and every NaN stays where it was**,
// which is neither what mean() nor what max() does and is measured against
// R2026a 2026-09-02: normalize([1 NaN 3]) is [-0.7071 NaN 0.7071], the
// z-score of 1 and 3 alone. A CONSTANT slice answers NaN for ZScore and
// Scale (a real 0/0) but `lo` everywhere for Range, the same rule
// rescale() already follows.
enum class Normalization { ZScore, Center, Scale, Range, Norm };
static ICoreMatrix normalize(const Normalization& how, const ICoreMatrix& a, int dim,
const ICoreMatrix* parameter, std::string* whyNot = nullptr);
// The mov* family. The window is described by how far it reaches BACK and
// FORWARD from each element; `window()` derives that pair from MATLAB's own
// argument, which is either a length k or an explicit [kb kf]. Endpoints
// shrink (MATLAB's default 'shrink' rule), so the first and last answers
// are computed over a shorter window rather than dropped.
enum class Moving { Mean, Sum, Median, Max, Min, Std, Var };
static ICoreMatrix moving(const Moving& which, const ICoreMatrix& a,
const long long& back, const long long& forward,
int dim, bool normalizeByN,
std::string* whyNot = nullptr);
// k odd is centred; k even reaches one further back than forward
// (kb = floor(k/2), kf = ceil(k/2) - 1, which is also what MATLAB does with
// a fractional k); [kb kf] says it outright.
static bool window(const ICoreMatrix& k, long long& back, long long& forward,
const std::string& who, std::string* whyNot = nullptr);
};
};ICoreText.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreText.h
ICoreText#
ICoreText.h:37 · class · 20 declaration(s)
Core MATLAB's text surface (the core-MATLAB board, C11.3-C11.12): the conversions between numbers and text, and the operations on a char array.
class ICoreText {
public:
// ---- sprintf (C11.3) ---------------------------------------------------
// MATLAB's sprintf. The arguments arrive FLATTENED into one stream, which
// is MATLAB's own model: every argument contributes its elements in
// column-major order and `isChar` says, per element, whether it came from
// a char array. The two vectors are parallel and must be the same length.
// `anyArguments` distinguishes "sprintf(fmt)" from "sprintf(fmt, [])" --
// see below, they differ.
//
// The conversions are `%d %i %u %o %x %X %f %e %E %g %G %s %c %%`, with the
// `-+ 0#` flags, a width, a precision, `*` taking either from the stream,
// and C length modifiers (`l`, `h`, `ll`) accepted and ignored. Escapes in
// the FORMAT are sprintf's own, not the lexer's: `\n \t \r \\ \a \b \f \v`
// and `%%`.
//
// THREE rules about running out of arguments, all three measured:
// 1. no arguments at all -- output stops at the first conversion, so
// sprintf("a%db") is "a";
// 2. arguments given but the stream is EMPTY (sprintf("a%db", [])) --
// the format is emitted once and each conversion contributes nothing,
// so that is "ab", NOT "a";
// 3. the stream runs dry part-way -- the format CYCLES until it does,
// then stops at the unsatisfied conversion keeping the literal text
// before it: sprintf("%d-%d,", [1 2 3]) is "1-2,3-".
//
// And one about the wrong type, which is MATLAB's silent switch rather
// than an error: an INTEGER conversion handed a non-integral value falls
// back to `%e` with the same flags, so sprintf("%d", 1.5) is
// "1.500000e+00". `%s` handed an integral number prints its character
// (sprintf("%s", 65) is "A") and a non-integral one falls back the same
// way; `%s` handed characters consumes the whole run of them.
static std::string format(const std::string& fmt,
const std::vector<double>& values,
const std::vector<bool>& isChar,
bool anyArguments,
std::string* whyNot = nullptr);
// ---- number to text (C11.4, C11.6, C11.12) -----------------------------
// num2str(x) and string(x) -- the same text for a numeric argument.
//
// MATLAB picks the format from the DATA, and the arithmetic is exact
// rather than approximate, because the padding between the entries of a
// vector is a visible part of the answer:
// - every finite entry integral and max|x| below flintmax: "%<w>d" with
// w = (digits in the largest magnitude) + 2;
// - otherwise "%<w>.<p>g" with p = floor(log10(max|x|)) + 5 clamped to
// [5, 16] and w = p + 7;
// - either way, w gains 2 more when any entry is Inf or NaN, and the
// finished row is trimmed at both ends.
// So num2str(pi) is "3.1416", num2str(123.456) is "123.456" (the width
// follows the magnitude, it is not a fixed 5 significant digits), and
// num2str([1 2.5]) is "1" then NINE spaces then "2.5".
//
// A char argument is answered unchanged, as MATLAB does. A matrix with
// more than one ROW is refused: MATLAB answers a char MATRIX there and
// this console has no such value (C0.3).
static std::string numberToText(const ICoreMatrix& a, std::string* whyNot = nullptr);
// num2str(x, n) -- n significant digits, "%<n+7>.<n>g" per entry, same
// trim. num2str(pi, 8) is "3.1415927".
static std::string numberToText(const ICoreMatrix& a, int significantDigits,
std::string* whyNot = nullptr);
// int2str(x) -- every entry rounded to the nearest integer, HALF AWAY FROM
// ZERO (int2str(2.5) is 3 and int2str(-2.5) is -3), then num2str's integer
// path. Inf and NaN pass through as words.
static std::string integerToText(const ICoreMatrix& a, std::string* whyNot = nullptr);
// mat2str(A) and mat2str(A, n) -- the form that reads BACK as an
// expression: "[1 2;3 4]", a scalar bare, an empty "[]". Default precision
// is 15 significant digits, which is what makes mat2str(pi) round-trip
// where num2str(pi) does not. A LOGICAL answers the words -- "true",
// "[true false]" -- because that is what reads back as a logical (C9.4).
static std::string matrixToText(const ICoreMatrix& a, int significantDigits);
// ---- text to number (C11.5) --------------------------------------------
// str2double(s). NaN for anything that is not ONE number, which is
// MATLAB's answer for "abc", for "" and for "1 2" alike. Surrounding
// whitespace is allowed, embedded commas are DROPPED (str2double("1,234.5")
// is 1234.5 and str2double("3,14") is 314, not 3.14), "Inf"/"NaN" are
// taken in any case with an optional sign, and a "0x" prefix is read as
// hexadecimal. There is no binary form: str2double("0b101") is NaN.
static double textToNumber(const std::string& s);
// ---- char array operations (C11.7-C11.12) ------------------------------
// strcat(...). MATLAB REMOVES THE TRAILING WHITESPACE of every char
// argument and keeps the leading whitespace, so strcat("a ", "b") is "ab"
// while strcat(" a", "b") is " ab". `[s1 s2]` does not do this -- that
// is plain concatenation and it is the reason both spellings exist.
static std::string concatenateTrimmed(const std::vector<std::string>& parts);
// strcmp / strcmpi / strncmp / strncmpi (C11.8). `n` compares the first n
// characters only and is false when EITHER side is shorter than n; n = 0
// is true.
enum class Case { Sensitive, Insensitive };
static bool same(const std::string& a, const std::string& b, const Case& how);
static bool samePrefix(const std::string& a, const std::string& b, size_t n,
const Case& how);
// upper / lower / strtrim / deblank / blanks (C11.9). strtrim takes
// whitespace off BOTH ends, deblank only off the end -- that is the whole
// difference between them.
static std::string upperCase(const std::string& s);
static std::string lowerCase(const std::string& s);
static std::string trimmed(const std::string& s);
static std::string trimmedEnd(const std::string& s);
// strfind(s, pattern) -- every 1-based start position, OVERLAPPING:
// strfind("aaa", "aa") is [1 2], not [1]. An empty pattern matches
// nothing.
static std::vector<double> findAll(const std::string& s, const std::string& pattern);
// strrep and replace (C11.10) -- and they DISAGREE, which is not a bug in
// either. strrep replaces at every position strfind reports, so the
// overlapping matches of strrep("aaaa", "aa", "b") each fire and the
// answer is "bbb"; replace scans left to right and skips past what it
// matched, so replace("aaaa", "aa", "b") is "bb". Both were measured.
static std::string replaceOverlapping(const std::string& s, const std::string& pattern,
const std::string& with);
static std::string replaceScanning(const std::string& s, const std::string& pattern,
const std::string& with);
// contains / startsWith / endsWith. An empty pattern is true for all three.
static bool contains(const std::string& s, const std::string& pattern);
static bool startsWith(const std::string& s, const std::string& pattern);
static bool endsWith(const std::string& s, const std::string& pattern);
// extractBefore / extractAfter, in both MATLAB spellings: before or after
// the FIRST occurrence of a pattern, or before or after a 1-based
// POSITION. A pattern that is not there answers "" for both, and so does a
// position at the far end.
static std::string before(const std::string& s, const std::string& pattern);
static std::string after(const std::string& s, const std::string& pattern);
static std::string beforePosition(const std::string& s, size_t position);
static std::string afterPosition(const std::string& s, size_t position);
// ---- regular expressions (C11.13) --------------------------------------
//
// The flavour is ECMAScript (std::regex), not MATLAB's own. The common
// subset -- character classes, quantifiers, anchors, groups, alternation,
// `\d \w \s`, `$1` backreferences in a replacement -- means the same thing
// in both, and the parity corpus only uses that subset. MATLAB-only
// constructs (named tokens `(?<name>...)`, lookbehind, `\<`) are a stated
// divergence rather than a silent one: a pattern this engine cannot read
// is REFUSED with the engine's own reason, never quietly matched
// differently.
//
// `$0` is translated to ECMAScript's `$&` on the way in, because MATLAB
// spells the whole match `$0` and ECMAScript does not have it at all --
// measured: regexprep("ab", "b", "$0$0") is "abb".
//
// matchPositions regexp(s, pat) / regexp(s, pat, "start") / "end" -- the
// 1-based positions of every match, as a row
// firstMatch regexp(s, pat, "once", "match") -- the matched TEXT, or
// "" when there is none (`found` says which)
// replaceAll regexprep(s, pat, rep)
static std::vector<double> matchPositions(const std::string& s, const std::string& pattern,
bool ignoreCase, bool wantEnd,
std::string* whyNot = nullptr);
static std::string firstMatch(const std::string& s, const std::string& pattern,
bool ignoreCase, bool* found, std::string* whyNot = nullptr);
static std::string replaceAll(const std::string& s, const std::string& pattern,
const std::string& with, bool ignoreCase,
std::string* whyNot = nullptr);
// char(A) (C11.12) -- code points to characters. The conversion TRUNCATES
// toward zero rather than rounding: char(65.9) is "A", not "B".
static std::string fromCodePoints(const ICoreMatrix& a, std::string* whyNot = nullptr);
// double(s) -- the code points of s as a 1 x n ROW, MATLAB's shape for a
// char array.
static ICoreMatrix toCodePoints(const std::string& s);
};
};ICoreTimeSeries.h#
src/ICoreBlocks/ICoreMath/Foundation/ICoreTimeSeries.h
ICoreTimeSeries#
ICoreTimeSeries.h:24 · class · pImpl · 20 declaration(s)
A sampled signal: N time stamps and an N x M block of values, one column per channel.
class ICoreTimeSeries {
public:
// NOT `= default` any more: the residue at the bottom of this class is a
// unique_ptr, and a defaulted default constructor would leave it NULL --
// which compiles, links, and dereferences null on first use.
ICoreTimeSeries();
// The series is stored in the variables space and handed around by value,
// so it stays copyable; the unique_ptr residue deletes the implicit copy
// operations, so all four are written out in the .cpp.
ICoreTimeSeries(const ICoreTimeSeries& other);
ICoreTimeSeries& operator=(const ICoreTimeSeries& other);
ICoreTimeSeries(ICoreTimeSeries&& other) noexcept;
ICoreTimeSeries& operator=(ICoreTimeSeries&& other) noexcept;
~ICoreTimeSeries();
// Single-channel construction from parallel vectors.
ICoreTimeSeries(std::vector<double> times, std::vector<double> values);
// Multi-channel construction. valuesRowMajor holds times.size() * channels
// entries, sample-major. Produces an invalid series (isValid() == false)
// rather than throwing if the counts disagree.
ICoreTimeSeries(std::vector<double> times, std::vector<double> valuesRowMajor, size_t channels);
// The console constructor: an N-element time vector (either orientation) and
// an N x M value matrix. A 1 x N value row vector is accepted as N samples of
// one channel, since that is how a time-shaped vector is usually typed.
// Returns false with `error` set on any size mismatch, never a partial series.
static bool fromMatrices(const ICoreMatrix& time, const ICoreMatrix& values,
ICoreTimeSeries& out, std::string& error);
void clear();
void append(const double& time, const double& value); // single channel
void appendSample(const double& time, const std::vector<double>& channelValues);
void reserve(const size_t& expectedSamples);
// Drops the oldest samples until at most maxSamples remain. A recorder on a
// long or infinite run would otherwise grow without bound.
void truncateToMostRecent(const size_t& maxSamples);
[[nodiscard]] bool isValid() const; // one value row per time stamp
[[nodiscard]] bool isEmpty() const;
[[nodiscard]] size_t sampleCount() const;
[[nodiscard]] size_t channelCount() const;
[[nodiscard]] const std::vector<double>& getTimes() const;
// Flat, sample-major. For a single-channel series this is simply the values
// in order; for several channels the caller must stride by channelCount().
[[nodiscard]] const std::vector<double>& getValuesRowMajor() const;
// One channel's samples in order. Empty when the index is out of range.
[[nodiscard]] std::vector<double> getChannel(const size_t& channelIndex) const;
// Mean sample interval, or 0 when there are fewer than two samples. Only
// meaningful for a fixed-step run -- see the class note.
[[nodiscard]] double meanSamplingTime() const;
[[nodiscard]] ICoreMatrix getTimesAsColumnVector() const; // N x 1
[[nodiscard]] ICoreMatrix getValuesAsMatrix() const; // N x M
// [time, ch0, ch1, ...] as an N x (1 + M) matrix -- exactly the shape
// ICoreChart::plot() reads as "column 0 is X, every later column is a line".
[[nodiscard]] ICoreMatrix getAsTimeValueMatrix() const;
private:
class Impl; // the two-line residue; state lives here
std::unique_ptr<Impl> impl;
};