codac 2.0.0
Loading...
Searching...
No Matches
codac2 Namespace Reference

Classes

struct  Color
 Color structure, in RGBA or HSVA format. More...
 
struct  ColorMap
 Represents a set of RGB or HSV values. More...
 
class  ConvexPolygon
 Represents a convex polygon defined by vertices enclosed in IntervalVectors. More...
 
class  CtcDist
 Implements the distance constraint on \(\mathbf{a}\in\mathbb{R}^2\), \(\mathbf{b}\in\mathbb{R}^2\) and \(d\in\mathbb{R}\) such that: More...
 
class  CtcGaussElim
 Contractor for a linear system of interval equations, based on the classical Gauss elimination procedure. More...
 
class  CtcGaussSeidel
 Contractor for a linear system of interval equations, using a fixed-point approach based on the Gauss Seidel method. More...
 
class  CtcLinearPrecond
 Contractor for a linear system of interval equations, using a preconditioning method before calling some provided contractor. More...
 
class  CtcPolar
 Implements the polar constraint on \(x\in\mathbb{R}\), \(y\in\mathbb{R}\), \(\rho\in\mathbb{R}^+\) and \(\theta\in\mathbb{R}\) such that: More...
 
class  DefaultFigure
 Default view class, used to manage the default figure. More...
 
class  Ellipsoid
 Ellipsoid representation. More...
 
class  ExprBase
 Abstract base class for representing an expression. More...
 
class  ExprID
 A class representing a unique identifier for expressions. More...
 
class  Figure2D
 Figure2D class, used for 2D display. More...
 
class  Figure2D_IPE
 IPE output class. More...
 
class  Figure2D_VIBes
 VIBes output class. More...
 
class  Figure2DInterface
 Interface for 2D figures. More...
 
class  Figure3D
 Figure3D class, used for 3D figures. More...
 
class  FunctionArgsList
 A container class to manage a collection of function arguments. More...
 
class  FunctionBase
 A base class for functions (either analytic functions, or set functions). More...
 
class  Interval
 Interval class, for representing closed and connected subsets of \(\mathbb{R}\). More...
 
class  IntvFullPivLU
 Full pivot LU decomposition for matrix of interval, based on Eigen decomposition. The decomposition is of the form M = P{-1} [L][U] Q{-1} where P and Q are permutation matrices, and [L] and [U] are lower and upper interval matrices ([L] diagonal is 1). More...
 
class  OperationExprBase
 A base class for expressions representing operations with multiple operands. More...
 
class  OutputFigure2D
 Output figure class, used to manage the output figures. More...
 
class  Polygon
 Represents a polygon (convex or non-convex) defined by its vertices enclosed in IntervalVectors. More...
 
class  Segment
 Represents a geometric segment defined by two points enclosed in IntervalVectors. More...
 
class  SepChi
 A separator on the Chi constraint If [x] in Sa, then Sb, else Sc. More...
 
struct  StyleProperties
 Style properties structure, to specify the style of a shape. More...
 
class  VarBase
 Abstract base class for representing variables in analytic or set functions. More...
 

Typedefs

using IntervalMatrix = Eigen::Matrix<Interval,-1,-1>
 Alias for a dynamic-size matrix of intervals.
 
using IntervalRow = Eigen::Matrix<Interval,1,-1>
 Alias for a dynamic-size row vector of intervals.
 
using IntervalVector = Eigen::Matrix<Interval,-1,1>
 Alias for a dynamic-size column vector of intervals.
 
using Matrix = Eigen::Matrix<double,-1,-1>
 Alias for a dynamic-size matrix of doubles.
 
using Row = Eigen::Matrix<double,1,-1>
 Alias for a dynamically-sized row vector of doubles.
 
using Vector = Eigen::Matrix<double,-1,1>
 Alias for a dynamically-sized column vector of doubles.
 

Enumerations

enum class  BoolInterval { }
 Enumeration representing a boolean interval. More...
 
enum class  OrientationInterval { }
 Enumeration representing feasible orientations. Can be used to assess an oriented angle, or the alignement of three points. More...
 
enum  Model
 Color model (RGB or HSV) More...
 

Functions

BoolInterval operator~ (BoolInterval x)
 Returns the complementary of a BoolInterval.
 
std::ostream & operator<< (std::ostream &os, const BoolInterval &x)
 Streams out a BoolInterval.
 
Ellipsoid operator+ (const Ellipsoid &e1, const Ellipsoid &e2)
 Compute the Minkowski sum of two ellipsoids.
 
Ellipsoid linear_mapping (const Ellipsoid &e, const Matrix &A, const Vector &b)
 Guaranteed linear evaluation A*e+b, considering the rounding errors.
 
Ellipsoid unreliable_linear_mapping (const Ellipsoid &e, const Matrix &A, const Vector &b)
 Nonrigorous linear evaluation A*e+b.
 
Ellipsoid nonlinear_mapping (const Ellipsoid &e, const AnalyticFunction< VectorType > &f)
 (Rigorous?) non-linear evaluation f(e)
 
Ellipsoid nonlinear_mapping (const Ellipsoid &e, const AnalyticFunction< VectorType > &f, const Vector &trig, const Vector &q)
 (Rigorous?) non-linear evaluation f(e), with parameters
 
Matrix nonlinear_mapping_base (const Matrix &G, const Matrix &J, const IntervalMatrix &J_box, const Vector &trig, const Vector &q)
 (Rigorous?) non-linear evaluation f(e), from Jacobian information
 
std::ostream & operator<< (std::ostream &os, const Ellipsoid &e)
 Streams out an Ellipsoid.
 
BoolInterval stability_analysis (const AnalyticFunction< VectorType > &f, unsigned int alpha_max, Ellipsoid &e, Ellipsoid &e_out, bool verbose=false)
 ...
 
Matrix solve_discrete_lyapunov (const Matrix &A, const Matrix &Q)
 ...
 
std::ostream & operator<< (std::ostream &os, const Interval &x)
 Streams out this.
 
Interval operator""_i (long double x)
 Codac defined literals allowing to produce an interval from a double.
 
double previous_float (double x)
 Returns the previous representable double-precision floating-point value before x.
 
double next_float (double x)
 Returns the next representable double-precision floating-point value after x.
 
Interval operator& (const Interval &x, const Interval &y)
 Returns the intersection of two intervals: \([x]\cap[y]\).
 
Interval operator| (const Interval &x, double y)
 Returns the squared-union of an interval and a real: \([x]\sqcup\{y\}\).
 
Interval operator| (const Interval &x, const Interval &y)
 Returns the squared-union of two intervals: \([x]\sqcup[y]\).
 
const Intervaloperator+ (const Interval &x)
 Returns this.
 
Interval operator+ (const Interval &x, double y)
 Returns \([x]+y\) with \(y\in\mathbb{R}\).
 
Interval operator+ (double x, const Interval &y)
 Returns \(x+[y]\) with \(x\in\mathbb{R}\).
 
Interval operator+ (const Interval &x, const Interval &y)
 Returns \([x]+[y]\).
 
Interval operator- (const Interval &x, double y)
 Returns \([x]-y\) with \(y\in\mathbb{R}\).
 
Interval operator- (double x, const Interval &y)
 Returns \(x-[y]\) with \(x\in\mathbb{R}\).
 
Interval operator- (const Interval &x, const Interval &y)
 Returns \([x]-[y]\).
 
Interval operator* (const Interval &x, double y)
 Returns \([x]*y\) with \(y\in\mathbb{R}\).
 
Interval operator* (double x, const Interval &y)
 Returns \(x*[y]\) with \(x\in\mathbb{R}\).
 
Interval operator* (const Interval &x, const Interval &y)
 Returns \([x]*[y]\).
 
Interval operator/ (const Interval &x, double y)
 Returns \([x]/y\) with \(y\in\mathbb{R}\).
 
Interval operator/ (double x, const Interval &y)
 Returns \(x/[y]\) with \(x\in\mathbb{R}\).
 
Interval operator/ (const Interval &x, const Interval &y)
 Returns \([x]/[y]\).
 
Interval sqr (const Interval &x)
 Returns \([x]^2\).
 
Interval sqrt (const Interval &x)
 Returns \(\sqrt{[x]}\).
 
Interval pow (const Interval &x, int n)
 Returns \([x]^n\), \(n\in\mathbb{Z}\).
 
Interval pow (const Interval &x, double p)
 Returns \([x]^p\), \(p\in\mathbb{R}\).
 
Interval pow (const Interval &x, const Interval &p)
 Returns \([x]^{[p]}\), \(p\in\mathbb{IR}\).
 
Interval root (const Interval &x, int p)
 Returns the p-th root: \(\sqrt[p]{[x]}\).
 
Interval exp (const Interval &x)
 Returns \(\exp([x])\).
 
Interval log (const Interval &x)
 Returns \(\log([x])\).
 
Interval cos (const Interval &x)
 Returns \(\cos([x])\).
 
Interval sin (const Interval &x)
 Returns \(\sin([x])\).
 
Interval tan (const Interval &x)
 Returns \(\tan([x])\).
 
Interval acos (const Interval &x)
 Returns \(\acos([x])\).
 
Interval asin (const Interval &x)
 Returns \(\asin([x])\).
 
Interval atan (const Interval &x)
 Returns \(\atan([x])\).
 
Interval atan2 (const Interval &y, const Interval &x)
 Returns \(\mathrm{arctan2}([y],[x])\).
 
Interval cosh (const Interval &x)
 Returns \(\cosh([x])\).
 
Interval sinh (const Interval &x)
 Returns \(\sinh([x])\).
 
Interval tanh (const Interval &x)
 Returns \(\tanh([x])\).
 
Interval acosh (const Interval &x)
 Returns \(\acosh([x])\).
 
Interval asinh (const Interval &x)
 Returns \(\asinh([x])\).
 
Interval atanh (const Interval &x)
 Returns \(\atanh([x])\).
 
Interval abs (const Interval &x)
 Returns \(\mid[x]\mid = \left\{\mid x \mid, x\in[x]\right\}\).
 
Interval min (const Interval &x, const Interval &y)
 Returns \(\min([x],[y])=\left\{\min(x,y), x\in[x], y\in[y]\right\}\).
 
Interval max (const Interval &x, const Interval &y)
 Returns \(\max([x],[y])=\left\{\max(x,y), x\in[x], y\in[y]\right\}\).
 
Interval sign (const Interval &x)
 Returns \(\sign([x])=\left[\left\{\sign(x), x\in[x]\right\}\right]\).
 
Interval integer (const Interval &x)
 Returns the largest integer interval included in \([x]\).
 
Interval floor (const Interval &x)
 Returns floor of \([x]\).
 
Interval ceil (const Interval &x)
 Returns ceil of \([x]\).
 
Interval chi (const Interval &x, const Interval &y, const Interval &z)
 Return \([y]\) if \(x^+ <= 0\), \([z]\) if \(x^- > 0\), \([y]\sqcup[z]\) else.
 
std::ostream & operator<< (std::ostream &os, const IntervalMatrix &x)
 Stream output operator for IntervalMatrix objects.
 
std::ostream & operator<< (std::ostream &os, const IntervalRow &x)
 Stream output operator for IntervalRow objects.
 
std::ostream & operator<< (std::ostream &os, const IntervalVector &x)
 Stream output operator for IntervalVector objects.
 
ConvexPolygon operator& (const ConvexPolygon &p1, const ConvexPolygon &p2)
 Computes the intersection of two convex polygons.
 
std::ostream & operator<< (std::ostream &os, const OrientationInterval &x)
 Streams out a OrientationInterval.
 
OrientationInterval orientation (const IntervalVector &p1, const IntervalVector &p2, const IntervalVector &p3)
 Computes the orientation of an ordered triplet of 2D points.
 
BoolInterval aligned (const IntervalVector &p1, const IntervalVector &p2, const IntervalVector &p3)
 Checks whether three 2D points are aligned (colinear).
 
std::vector< IntervalVectorconvex_hull (std::vector< IntervalVector > pts)
 Computes the convex hull of a set of 2d points.
 
IntervalVector rotate (const IntervalVector &p, const Interval &a, const IntervalVector &c=Vector::zero(2))
 Rotates a 2D interval vector around a point by a given angle.
 
template<typename T>
rotate (const T &x, const Interval &a, const IntervalVector &c=Vector::zero(2))
 Rotates a 2D object (Segment, Polygon, etc) around a point by a given angle.
 
std::ostream & operator<< (std::ostream &str, const Polygon &p)
 Stream output operator for Polygon.
 
IntervalVector operator& (const Segment &e1, const Segment &e2)
 Computes the intersection of two segments.
 
IntervalVector proj_intersection (const Segment &e1, const Segment &e2)
 Computes the projected intersection of two segments.
 
BoolInterval colinear (const Segment &e1, const Segment &e2)
 Checks if two segments are colinear.
 
std::ostream & operator<< (std::ostream &str, const Segment &e)
 Stream output operator for Segment.
 
template<LeftOrRightInv O = LEFT_INV, typename OtherDerived, typename OtherDerived_>
IntervalMatrix inverse_correction (const Eigen::MatrixBase< OtherDerived > &A, const Eigen::MatrixBase< OtherDerived_ > &B)
 Correct the approximation of the inverse \(\mathbf{B}\approx\mathbf{A}^{-1}\) of a square matrix \(\mathbf{A}\) by providing a reliable enclosure \([\mathbf{A}^{-1}]\).
 
template<typename OtherDerived>
IntervalMatrix inverse_enclosure (const Eigen::MatrixBase< OtherDerived > &A)
 Enclosure of the inverse of a (non-singular) matrix expression, possibly an interval matrix.
 
IntervalMatrix infinite_sum_enclosure (const IntervalMatrix &A, double &mrad)
 Compute an upper bound of \(\left([\mathbf{A}]+[\mathbf{A}]^2+[\mathbf{A}]^3+\dots\right)\), with \([\mathbf{A}]\) a matrix of intervals as an "error term" (uses only bounds on coefficients).
 
template<typename OtherDerived>
auto abs (const Eigen::MatrixBase< OtherDerived > &x)
 Compute the element-wise absolute value of a matrix.
 
template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto floor (const Eigen::MatrixBase< OtherDerived > &x)
 Compute the element-wise floor of a matrix.
 
template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto ceil (const Eigen::MatrixBase< OtherDerived > &x)
 Compute the element-wise ceiling of a matrix.
 
template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto round (const Eigen::MatrixBase< OtherDerived > &x)
 Compute the element-wise rounding of a matrix.
 
Eigen::IOFormat codac_row_fmt ()
 Provides an Eigen IOFormat for formatting row vectors.
 
Eigen::IOFormat codac_vector_fmt ()
 Provides an Eigen IOFormat for formatting column vectors.
 
Eigen::IOFormat codac_matrix_fmt ()
 Provides an Eigen IOFormat for formatting matrices.
 
std::ostream & operator<< (std::ostream &os, const Matrix &x)
 Stream output operator for Matrix objects.
 
std::ostream & operator<< (std::ostream &os, const Row &x)
 Stream output operator for Row objects.
 
std::ostream & operator<< (std::ostream &os, const Vector &x)
 Stream output operator for Vector objects.
 
template<typename T>
requires (std::is_trivially_copyable_v<T>)
void serialize (std::ostream &f, const T &x)
 Writes the binary representation of a trivially copyable object to the given output stream.
 
template<typename T>
requires (std::is_trivially_copyable_v<T>)
void deserialize (std::istream &f, T &x)
 Reads the binary representation of a trivially copyable object from the given input stream.
 
void serialize (std::ostream &f, const Interval &x)
 Writes the binary representation of an Interval object to the given output stream.
 
void deserialize (std::istream &f, Interval &x)
 Creates an Interval object from the binary representation given in an input stream.
 
template<typename T, int R = -1, int C = -1>
void serialize (std::ostream &f, const Eigen::Matrix< T, R, C > &x)
 Writes the binary representation of an Eigen::Matrix to the given output stream. The structure can be a matrix, a row, a vector, with double or Interval components.
 
template<typename T, int R = -1, int C = -1>
void deserialize (std::istream &f, Eigen::Matrix< T, R, C > &x)
 Reads the binary representation of an Eigen::Matrix from the given input stream. The structure can be a matrix, a row, a vector, with double or Interval components.
 
template<typename T>
void serialize (std::ostream &f, const SampledTraj< T > &x)
 Writes the binary representation of a SampledTraj object to the given output stream.
 
template<typename T>
void deserialize (std::istream &f, SampledTraj< T > &x)
 Reads the binary representation of a SampledTraj object from the given input stream.
 
template<typename T>
const SampledTraj< T > & operator+ (const SampledTraj< T > &x1)
 \(x1(\cdot)\)
 
template<typename T>
SampledTraj< T > operator+ (const SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_binary_traj_traj(operator_add< T >)
 \(x_1(\cdot)+x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > operator+ (const SampledTraj< T > &x1, const Q &x2) macro_binary_traj_real(operator_add< T >)
 \(x_1(\cdot)+x_2\)
 
template<typename T, typename Q>
SampledTraj< T > operator+ (const Q &x1, const SampledTraj< T > &x2) macro_binary_real_traj(operator_add< T >)
 \(x+x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > & operator+= (SampledTraj< T > &x1, const Q &x2) macro_member_binary_traj_real(operator_add< T >)
 Operates +=.
 
template<typename T>
SampledTraj< T > & operator+= (SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_member_binary_traj_traj(operator_add< T >)
 Operates +=.
 
template<typename T>
SampledTraj< T > operator- (const SampledTraj< T > &x1)
 \(-x_1(\cdot)\)
 
template<typename T>
SampledTraj< T > operator- (const SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_binary_traj_traj(operator_sub< T >)
 \(x_1(\cdot)-x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > operator- (const SampledTraj< T > &x1, const Q &x2) macro_binary_traj_real(operator_sub< T >)
 \(x_1(\cdot)-x_2\)
 
template<typename T, typename Q>
SampledTraj< T > operator- (const Q &x1, const SampledTraj< T > &x2) macro_binary_real_traj(operator_sub< T >)
 \(x-x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > & operator-= (SampledTraj< T > &x1, const Q &x2) macro_member_binary_traj_real(operator_sub< T >)
 Operates -=.
 
template<typename T>
SampledTraj< T > & operator-= (SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_member_binary_traj_traj(operator_sub< T >)
 Operates -=.
 
template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > operator* (double x1, const SampledTraj< T > &x2) macro_binary_real_traj(operator_mul_scal< T >)
 \(x_1\cdot x_2(\cdot)\)
 
template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > operator* (const SampledTraj< T > &x1, double x2) macro_binary_traj_real(operator_mul_scal< T >)
 \(x_1(\cdot)\cdot x_2\)
 
template<typename T>
SampledTraj< T > operator* (const SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_binary_traj_traj(operator_mul< T >)
 \(x_1(\cdot)\cdot x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > operator* (const SampledTraj< T > &x1, const Q &x2) macro_binary_traj_real(operator_mul< T >)
 \(x_1(\cdot)\cdot x_2\)
 
template<typename T, typename Q>
SampledTraj< T > operator* (const Q &x1, const SampledTraj< T > &x2) macro_binary_real_traj(operator_mul< T >)
 \(x\cdot x_2(\cdot)\)
 
SampledTraj< Vectoroperator* (const SampledTraj< Matrix > &x1, const SampledTraj< Vector > &x2) macro_binary_traj_traj(operator_mul_vec)
 \(x_1(\cdot)\cdot x_2\)
 
template<typename T, typename Q>
SampledTraj< T > & operator*= (SampledTraj< T > &x1, const Q &x2) macro_member_binary_traj_real(operator_mul< T >)
 Operates *=.
 
template<typename T>
SampledTraj< T > & operator*= (SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_member_binary_traj_traj(operator_mul< T >)
 Operates *=.
 
template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > operator/ (const SampledTraj< T > &x1, double x2) macro_binary_traj_real(operator_div_scal< T >)
 \(x_2(\cdot)/x_1\)
 
template<typename T>
SampledTraj< T > operator/ (const SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_binary_traj_traj(operator_div< T >)
 \(x_1(\cdot)/x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > operator/ (const SampledTraj< T > &x1, const Q &x2) macro_binary_traj_real(operator_div< T >)
 \(x_1(\cdot)/x_2\)
 
template<typename T, typename Q>
SampledTraj< T > operator/ (const Q &x1, const SampledTraj< T > &x2) macro_binary_real_traj(operator_div< T >)
 \(x/x_2(\cdot)\)
 
template<typename T, typename Q>
SampledTraj< T > & operator/= (SampledTraj< T > &x1, const Q &x2) macro_member_binary_traj_real(operator_div< T >)
 Operates /=.
 
template<typename T>
SampledTraj< T > & operator/= (SampledTraj< T > &x1, const SampledTraj< T > &x2) macro_member_binary_traj_traj(operator_div< T >)
 Operates /=.
 
SampledTraj< double > sqr (const SampledTraj< double > &x1)
 \(x^2(\cdot)\)
 
SampledTraj< double > sqrt (const SampledTraj< double > &x1)
 \(\sqrt{x_1(\cdot)}\)
 
SampledTraj< double > pow (const SampledTraj< double > &x1, int x2)
 \(x^x_2(\cdot)\)
 
SampledTraj< double > pow (const SampledTraj< double > &x1, double x2)
 \(x^x_2(\cdot)\)
 
SampledTraj< double > root (const SampledTraj< double > &x1, int x2)
 \(x^x_2(\cdot)\)
 
SampledTraj< double > exp (const SampledTraj< double > &x1)
 \(\exp(x_1(\cdot))\)
 
SampledTraj< double > log (const SampledTraj< double > &x1)
 \(\log(x_1(\cdot))\)
 
SampledTraj< double > cos (const SampledTraj< double > &x1)
 \(\cos(x_1(\cdot))\)
 
SampledTraj< double > sin (const SampledTraj< double > &x1)
 \(\sin(x_1(\cdot))\)
 
SampledTraj< double > tan (const SampledTraj< double > &x1)
 \(\tan(x_1(\cdot))\)
 
SampledTraj< double > acos (const SampledTraj< double > &x1)
 \(\arccos(x_1(\cdot))\)
 
SampledTraj< double > asin (const SampledTraj< double > &x1)
 \(\arcsin(x_1(\cdot))\)
 
SampledTraj< double > atan (const SampledTraj< double > &x1)
 \(\arctan(x_1(\cdot))\)
 
SampledTraj< double > atan2 (const SampledTraj< double > &x1, const SampledTraj< double > &x2)
 \(\mathrm{arctan2}(x_1(\cdot),x_2(\cdot))\)
 
SampledTraj< double > atan2 (const SampledTraj< double > &x1, double x2)
 \(\mathrm{arctan2}(x_1(\cdot),x_2)\)
 
SampledTraj< double > atan2 (double x1, const SampledTraj< double > &x2)
 \(\mathrm{arctan2}(x_1, x_2(\cdot))\)
 
SampledTraj< double > cosh (const SampledTraj< double > &x1)
 \(\cosh(x_1(\cdot))\)
 
SampledTraj< double > sinh (const SampledTraj< double > &x1)
 \(\sinh(x_1(\cdot))\)
 
SampledTraj< double > tanh (const SampledTraj< double > &x1)
 \(\tanh(x_1(\cdot))\)
 
SampledTraj< double > acosh (const SampledTraj< double > &x1)
 \(\mathrm{arccosh}(x_1(\cdot))\)
 
SampledTraj< double > asinh (const SampledTraj< double > &x1)
 \(\mathrm{arcsinh}(x_1(\cdot))\)
 
SampledTraj< double > atanh (const SampledTraj< double > &x1)
 \(\mathrm{arctanh}(x_1(\cdot))\)
 
SampledTraj< double > abs (const SampledTraj< double > &x1)
 \(\mid x_1(\cdot)\mid\)
 
SampledTraj< double > min (const SampledTraj< double > &x1, const SampledTraj< double > &x2)
 \(\min(x_1(\cdot),x_2(\cdot))\)
 
SampledTraj< double > min (const SampledTraj< double > &x1, double x2)
 \(\min(x_1(\cdot),x_2)\)
 
SampledTraj< double > min (double x1, const SampledTraj< double > &x2)
 \(\min(x_1, x_2(\cdot))\)
 
SampledTraj< double > max (const SampledTraj< double > &x1, const SampledTraj< double > &x2)
 \(\max(x_1(\cdot),x_2(\cdot))\)
 
SampledTraj< double > max (const SampledTraj< double > &x1, double x2)
 \(\max(x_1(\cdot),x_2)\)
 
SampledTraj< double > max (double x1, const SampledTraj< double > &x2)
 \(\max(x_1, x_2(\cdot))\)
 
capd::Interval to_capd (const codac2::Interval &x)
 Converts a Codac Interval object into an CAPD Interval object.
 
codac2::Interval to_codac (const capd::Interval &x)
 Converts an CAPD Interval object into a Codac Interval object.
 
capd::IVector to_capd (const codac2::IntervalVector &x)
 Converts a Codac IntervalVector object into an CAPD IVector object.
 
codac2::IntervalVector to_codac (const capd::IVector &x)
 Converts an CAPD IVector object into a Codac IntervalVector object.
 
capd::IMatrix to_capd (const codac2::IntervalMatrix &x)
 Converts a Codac IntervalMatrix object into an CAPD IMatrix object.
 
codac2::IntervalMatrix to_codac (const capd::IMatrix &x)
 Converts an CAPD IMatrix object into a Codac IntervalMatrix object.
 
void draw_while_paving (const IntervalVector &x0, std::shared_ptr< const CtcBase< IntervalVector > > c, double eps, std::shared_ptr< Figure2D > fig=nullptr)
 Draws the paving while it is being computed.
 
void draw_while_paving (const IntervalVector &x0, const CtcBase< IntervalVector > &c, double eps, std::shared_ptr< Figure2D > fig=nullptr)
 Draws the paving while it is being computed.
 
void draw_while_paving (const IntervalVector &x0, std::shared_ptr< const SepBase > s, double eps, std::shared_ptr< Figure2D > fig=nullptr)
 Draws the paving while it is being computed.
 
void draw_while_paving (const IntervalVector &x0, const SepBase &s, double eps, std::shared_ptr< Figure2D > fig=nullptr)
 Draws the paving while it is being computed.
 

Detailed Description

codac2_ibex_impl.h

Date
2024
Author
Gilles Chabert, (Simon Rohou)
License: GNU Lesser General Public License (LGPL)

Typedef Documentation

◆ IntervalMatrix

using codac2::IntervalMatrix = Eigen::Matrix<Interval,-1,-1>

Alias for a dynamic-size matrix of intervals.

Represents a matrix with a dynamic number of rows and columns, where each element is an Interval object.

This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<Interval,-1,-1>.

◆ IntervalRow

using codac2::IntervalRow = Eigen::Matrix<Interval,1,-1>

Alias for a dynamic-size row vector of intervals.

Represents a row vector with a dynamic number of columns, where each element is an Interval object.

This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<Interval,1,-1>.

◆ IntervalVector

using codac2::IntervalVector = Eigen::Matrix<Interval,-1,1>

Alias for a dynamic-size column vector of intervals.

Represents a column vector with a dynamic number of rows, where each element is an Interval object.

This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<Interval,-1,1>.

◆ Matrix

using codac2::Matrix = Eigen::Matrix<double,-1,-1>

Alias for a dynamic-size matrix of doubles.

Represents a matrix with dynamic number of rows and columns, where each element is a double precision floating point number.

This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<double,-1,-1>.

See also
Vector, Row, Eigen::Matrix

◆ Row

using codac2::Row = Eigen::Matrix<double,1,-1>

Alias for a dynamically-sized row vector of doubles.

Defines a convenient shorthand for representing row vectors with dynamic size. This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<double,1,-1>.

See also
Eigen::Matrix

◆ Vector

using codac2::Vector = Eigen::Matrix<double,-1,1>

Alias for a dynamically-sized column vector of doubles.

Defines a convenient shorthand for representing column vectors with dynamic size. This type alias is based on Eigen's matrix template and corresponds to Eigen::Matrix<double,-1,1>.

See also
Eigen::Matrix

Enumeration Type Documentation

◆ BoolInterval

enum class codac2::BoolInterval
strong

Enumeration representing a boolean interval.

The logical operators & and | can be used to combine BoolInterval values.

Enumerator
EMPTY 

EMPTY is equivalent to the operation TRUE & FALSE.

UNKNOWN 

UNKNOWN is equivalent to the operation TRUE | FALSE.

23 {
24 FALSE = 0x01,
25 TRUE = 0x02,
27 EMPTY = 0x00,
29 UNKNOWN = 0x01 | 0x02
30 };
@ UNKNOWN
Definition codac2_BoolInterval.h:29
@ EMPTY
Definition codac2_BoolInterval.h:27

◆ OrientationInterval

enum class codac2::OrientationInterval
strong

Enumeration representing feasible orientations. Can be used to assess an oriented angle, or the alignement of three points.

The logical operator | can be used to combine OrientationInterval values.

Enumerator
UNKNOWN 

UNKNOWN: the orientation can be either clockwise, counter-clockwise, or aligned.

25 {
26 COLINEAR = 0x01,
27 CLOCKWISE = 0x02,
28 COUNTERCLOCKWISE = 0x04,
29 EMPTY = 0x00,
31 UNKNOWN = 0x01 | 0x02 | 0x04
32 };

◆ Model

Color model (RGB or HSV)

This enum is used to specify the color model of a color, either RGB or HSV.

28{ RGB, HSV };

Function Documentation

◆ operator~()

BoolInterval codac2::operator~ ( BoolInterval x)
inline

Returns the complementary of a BoolInterval.

Parameters
xthe boolean interval
Returns
the complementary
45 {
46 switch(x)
47 {
48 case BoolInterval::FALSE:
49 return BoolInterval::TRUE;
50 case BoolInterval::TRUE:
51 return BoolInterval::FALSE;
52 default:
54 }
55 }

◆ operator<<() [1/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const BoolInterval & x )
inline

Streams out a BoolInterval.

Parameters
osthe stream to be updated
xthe boolean interval to stream out
Returns
a reference to the updated stream
65 {
66 switch(x)
67 {
69 os << "[ empty ]";
70 break;
71 case BoolInterval::FALSE:
72 os << "[ false ]";
73 break;
74 case BoolInterval::TRUE:
75 os << "[ true ]";
76 break;
78 os << "[ true, false ]";
79 break;
80 }
81 return os;
82 }

◆ operator+() [1/9]

Ellipsoid codac2::operator+ ( const Ellipsoid & e1,
const Ellipsoid & e2 )

Compute the Minkowski sum of two ellipsoids.

Parameters
e1first ellipsoid
e2second ellipsoid
Returns
the Minkowski sum

◆ linear_mapping()

Ellipsoid codac2::linear_mapping ( const Ellipsoid & e,
const Matrix & A,
const Vector & b )

Guaranteed linear evaluation A*e+b, considering the rounding errors.

Parameters
einput ellipsoid
Amatrix
bvector
Returns
a rigorous outer enclosure of the linear mapping

◆ unreliable_linear_mapping()

Ellipsoid codac2::unreliable_linear_mapping ( const Ellipsoid & e,
const Matrix & A,
const Vector & b )

Nonrigorous linear evaluation A*e+b.

Note
This function is used in linear_mapping() and provides a faster output than its guaranteed counterpart.
Parameters
einput ellipsoid
Amatrix
bvector
Returns
a nonrigorous approximation of the linear mapping

◆ nonlinear_mapping() [1/2]

Ellipsoid codac2::nonlinear_mapping ( const Ellipsoid & e,
const AnalyticFunction< VectorType > & f )

(Rigorous?) non-linear evaluation f(e)

Parameters
einput ellipsoid
fnon-linear analytical function
Returns
a (rigorous?) outer enclosure of the non-linear mapping

◆ nonlinear_mapping() [2/2]

Ellipsoid codac2::nonlinear_mapping ( const Ellipsoid & e,
const AnalyticFunction< VectorType > & f,
const Vector & trig,
const Vector & q )

(Rigorous?) non-linear evaluation f(e), with parameters

Parameters
einput ellipsoid
fnon-linear analytical function
trig(?)
q(?)
Returns
a (rigorous?) outer enclosure of the non-linear mapping

◆ nonlinear_mapping_base()

Matrix codac2::nonlinear_mapping_base ( const Matrix & G,
const Matrix & J,
const IntervalMatrix & J_box,
const Vector & trig,
const Vector & q )

(Rigorous?) non-linear evaluation f(e), from Jacobian information

Parameters
G(?)
Japproximated Jacobian matrix of f
J_boxreliable enclosure of the Jacobian matrix of f
trig(?)
q(?)
Returns
a (rigorous?) outer enclosure of the non-linear mapping

◆ operator<<() [2/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const Ellipsoid & e )

Streams out an Ellipsoid.

Parameters
osthe stream to be updated
ethe ellipsoid stream out
Returns
a reference to the updated stream

◆ stability_analysis()

BoolInterval codac2::stability_analysis ( const AnalyticFunction< VectorType > & f,
unsigned int alpha_max,
Ellipsoid & e,
Ellipsoid & e_out,
bool verbose = false )

...

Parameters
f...
alpha_max...
e...
e_out...
verbose...
Returns
...

◆ solve_discrete_lyapunov()

Matrix codac2::solve_discrete_lyapunov ( const Matrix & A,
const Matrix & Q )

...

Parameters
A...
Q...
Returns
...

◆ operator<<() [3/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const Interval & x )
inline

Streams out this.

Parameters
osthe stream to be updated
xthe interval to stream out
Returns
a reference to the updated stream
589 {
590 gaol::interval::precision(os.precision());
591 gaol::operator<<(os,x);
592 return os;
593 }

◆ operator""_i()

Interval codac2::operator""_i ( long double x)
inline

Codac defined literals allowing to produce an interval from a double.

Parameters
xvalue to convert into an interval object
Returns
the interval \([x]\)
600 {
601 return Interval(x);
602 }
Interval class, for representing closed and connected subsets of .
Definition codac2_Interval.h:49

◆ previous_float()

double codac2::previous_float ( double x)
inline

Returns the previous representable double-precision floating-point value before x.

This function computes the largest double value that is strictly less than the input value x, effectively moving one step down in the floating-point representation.

See also
For obtaining the next representable double value greater than x, use next_float().
Parameters
xThe input double value.
Returns
The previous representable double value before x.
605 {
606 return gaol::previous_float(x);
607 }

◆ next_float()

double codac2::next_float ( double x)
inline

Returns the next representable double-precision floating-point value after x.

This function computes the smallest double value that is strictly greater than the input value x, effectively moving one step up in the floating-point representation.

See also
For obtaining the previous representable double value less than x, use previous_float().
Parameters
xThe input double value.
Returns
The next representable double value after x.
610 {
611 return gaol::next_float(x);
612 }

◆ operator&() [1/3]

Interval codac2::operator& ( const Interval & x,
const Interval & y )
inline

Returns the intersection of two intervals: \([x]\cap[y]\).

Note
Returns an empty interval if there is no intersection.
Parameters
xinterval value
yinterval value
Returns
intersection result
366 {
367 if(x.is_empty() || y.is_empty() || x.ub() < y.lb())
368 return Interval::empty();
369
370 else
371 return gaol::operator&(x,y);
372 }
bool is_empty() const
Tests if this is empty.
Definition codac2_Interval_impl.h:195
double ub() const
Returns the upper bound of this.
Definition codac2_Interval_impl.h:115
double lb() const
Returns the lower bound of this.
Definition codac2_Interval_impl.h:110
static Interval empty()
Provides an empty interval.
Definition codac2_Interval_impl.h:558

◆ operator|() [1/2]

Interval codac2::operator| ( const Interval & x,
double y )
inline

Returns the squared-union of an interval and a real: \([x]\sqcup\{y\}\).

Note
The squared-union is defined as: \([x]\sqcup[y]=\left[[x]\cup[y]\right]\)
Parameters
xinterval value
yreal value
Returns
squared-union result
375 {
376 return gaol::operator|(x,Interval(y));
377 }

◆ operator|() [2/2]

Interval codac2::operator| ( const Interval & x,
const Interval & y )
inline

Returns the squared-union of two intervals: \([x]\sqcup[y]\).

Note
The squared-union is defined as: \([x]\sqcup[y]=\left[[x]\cup[y]\right]\)
Parameters
xinterval value
yinterval value
Returns
squared-union result
380 {
381 return gaol::operator|(x,y);
382 }

◆ operator+() [2/9]

const Interval & codac2::operator+ ( const Interval & x)
inline

Returns this.

Note
This operator is only provided for consistency purposes.
Parameters
xinterval value
Returns
the same interval
385 {
386 return x;
387 }

◆ operator+() [3/9]

Interval codac2::operator+ ( const Interval & x,
double y )
inline

Returns \([x]+y\) with \(y\in\mathbb{R}\).

Parameters
xinterval value
yreal value
Returns
the addition result
390 {
391 if(y == -oo || y == oo)
392 return Interval::empty();
393
394 else
395 return gaol::operator+(x,y);
396 }

◆ operator+() [4/9]

Interval codac2::operator+ ( double x,
const Interval & y )
inline

Returns \(x+[y]\) with \(x\in\mathbb{R}\).

Parameters
xreal value
yinterval value
Returns
the addition result
399 {
400 if(x == -oo || x == oo)
401 return Interval::empty();
402
403 else
404 return gaol::operator+(x,y);
405 }

◆ operator+() [5/9]

Interval codac2::operator+ ( const Interval & x,
const Interval & y )
inline

Returns \([x]+[y]\).

Parameters
xinterval value
yinterval value
Returns
the addition result
408 {
409 return gaol::operator+(x,y);
410 }

◆ operator-() [1/7]

Interval codac2::operator- ( const Interval & x,
double y )
inline

Returns \([x]-y\) with \(y\in\mathbb{R}\).

Parameters
xinterval value
yreal value
Returns
the substraction result
413 {
414 if(y == -oo || y == oo)
415 return Interval::empty();
416
417 else
418 return gaol::operator-(x, y);
419 }

◆ operator-() [2/7]

Interval codac2::operator- ( double x,
const Interval & y )
inline

Returns \(x-[y]\) with \(x\in\mathbb{R}\).

Parameters
xreal value
yinterval value
Returns
the substraction result
422 {
423 if(x == -oo || x == oo)
424 return Interval::empty();
425
426 else
427 return gaol::operator-(x, y);
428 }

◆ operator-() [3/7]

Interval codac2::operator- ( const Interval & x,
const Interval & y )
inline

Returns \([x]-[y]\).

Parameters
xinterval value
yinterval value
Returns
the substraction result
431 {
432 return gaol::operator-(x, y);
433 }

◆ operator*() [1/9]

Interval codac2::operator* ( const Interval & x,
double y )
inline

Returns \([x]*y\) with \(y\in\mathbb{R}\).

Parameters
xinterval value
yreal value
Returns
the multiplication result
436 {
437 if(y == -oo || y == oo)
438 return Interval::empty();
439
440 else
441 return gaol::operator*(x,y);
442 }

◆ operator*() [2/9]

Interval codac2::operator* ( double x,
const Interval & y )
inline

Returns \(x*[y]\) with \(x\in\mathbb{R}\).

Parameters
xreal value
yinterval value
Returns
the multiplication result
445 {
446 if(x == -oo || x == oo)
447 return Interval::empty();
448
449 else
450 return gaol::operator*(x,y);
451 }

◆ operator*() [3/9]

Interval codac2::operator* ( const Interval & x,
const Interval & y )
inline

Returns \([x]*[y]\).

Parameters
xinterval value
yinterval value
Returns
the multiplication result
454 {
455 return gaol::operator*(x,y);
456 }

◆ operator/() [1/7]

Interval codac2::operator/ ( const Interval & x,
double y )
inline

Returns \([x]/y\) with \(y\in\mathbb{R}\).

Parameters
xinterval value
yreal value
Returns
the division result
459 {
460 if(y == -oo || y == oo)
461 return Interval::empty();
462
463 else
464 return gaol::operator/(x,y);
465 }

◆ operator/() [2/7]

Interval codac2::operator/ ( double x,
const Interval & y )
inline

Returns \(x/[y]\) with \(x\in\mathbb{R}\).

Parameters
xreal value
yinterval value
Returns
the division result
468 {
469 if(x == -oo || x == oo)
470 return Interval::empty();
471
472 else
473 return gaol::operator/(x,y);
474 }

◆ operator/() [3/7]

Interval codac2::operator/ ( const Interval & x,
const Interval & y )
inline

Returns \([x]/[y]\).

Parameters
xinterval value
yinterval value
Returns
the division result
477 {
478 return gaol::operator/(x,y);
479 }

◆ sqr() [1/2]

Interval codac2::sqr ( const Interval & x)
inline

Returns \([x]^2\).

Parameters
xinterval value
Returns
the operation result
22 {
23 return gaol::sqr(x);
24 }

◆ sqrt() [1/2]

Interval codac2::sqrt ( const Interval & x)
inline

Returns \(\sqrt{[x]}\).

Parameters
xinterval value
Returns
the operation result
27 {
28 Interval y = gaol::sqrt(x);
29 gaol::round_upward();
30 return y;
31 }

◆ pow() [1/5]

Interval codac2::pow ( const Interval & x,
int n )
inline

Returns \([x]^n\), \(n\in\mathbb{Z}\).

Parameters
xinterval value
ninteger power value
Returns
the operation result
34 {
35 Interval y = gaol::pow(x,p);
36 //gaol::round_upward(); // not necessary?
37 return y;
38 }

◆ pow() [2/5]

Interval codac2::pow ( const Interval & x,
double p )
inline

Returns \([x]^p\), \(p\in\mathbb{R}\).

Parameters
xinterval value
preal power value
Returns
the operation result
41 {
42 if(p == -oo || p == oo)
43 return Interval::empty();
44
45 else
46 {
47 Interval y = gaol::pow(x,p);
48 gaol::round_upward();
49 return y;
50 }
51 }

◆ pow() [3/5]

Interval codac2::pow ( const Interval & x,
const Interval & p )
inline

Returns \([x]^{[p]}\), \(p\in\mathbb{IR}\).

Parameters
xinterval value
pinterval power value
Returns
the operation result
54 {
55 Interval y = gaol::pow(x,p);
56 gaol::round_upward();
57 return y;
58 }

◆ root() [1/2]

Interval codac2::root ( const Interval & x,
int p )
inline

Returns the p-th root: \(\sqrt[p]{[x]}\).

Parameters
xinterval value
pinteger root
Returns
the operation result
61 {
62 // Get the root of the positive part (gaol does
63 // not consider negative values to be in the definition
64 // domain of the root function)
65
66 gaol::interval y = gaol::nth_root(x, p>=0 ? p : -p);
67
68 if(p%2 == 1 && x.lb() < 0)
69 y |= -gaol::nth_root(-x, p >= 0 ? p : -p);
70
71 if(p < 0)
72 y = 1.0/y;
73
74 gaol::round_upward();
75 return y;
76 }

◆ exp() [1/2]

Interval codac2::exp ( const Interval & x)
inline

Returns \(\exp([x])\).

Parameters
xinterval value
Returns
the operation result
79 {
80 Interval y = gaol::exp(x);
81 gaol::round_upward();
82 return y;
83 }

◆ log() [1/2]

Interval codac2::log ( const Interval & x)
inline

Returns \(\log([x])\).

Parameters
xinterval value
Returns
the operation result
86 {
87 if(x.ub() <= 0) // gaol returns (-oo,-DBL_MAX) if x.ub()==0, instead of empty set
88 return Interval::empty();
89
90 else
91 {
92 Interval y = gaol::log(x);
93 gaol::round_upward();
94 return y;
95 }
96 }

◆ cos() [1/2]

Interval codac2::cos ( const Interval & x)
inline

Returns \(\cos([x])\).

Parameters
xinterval value
Returns
the operation result
99 {
100 Interval y = gaol::cos(x);
101 gaol::round_upward();
102 return y;
103 }

◆ sin() [1/2]

Interval codac2::sin ( const Interval & x)
inline

Returns \(\sin([x])\).

Parameters
xinterval value
Returns
the operation result
106 {
107 Interval y = gaol::sin(x);
108 gaol::round_upward();
109 return y;
110 }

◆ tan() [1/2]

Interval codac2::tan ( const Interval & x)
inline

Returns \(\tan([x])\).

Parameters
xinterval value
Returns
the operation result
113 {
114 Interval y = gaol::tan(x);
115 gaol::round_upward();
116 return y;
117 }

◆ acos() [1/2]

Interval codac2::acos ( const Interval & x)
inline

Returns \(\acos([x])\).

Parameters
xinterval value
Returns
the operation result
120 {
121 Interval y = gaol::acos(x);
122 gaol::round_upward();
123 return y;
124 }

◆ asin() [1/2]

Interval codac2::asin ( const Interval & x)
inline

Returns \(\asin([x])\).

Parameters
xinterval value
Returns
the operation result
127 {
128 Interval y = gaol::asin(x);
129 gaol::round_upward();
130 return y;
131 }

◆ atan() [1/2]

Interval codac2::atan ( const Interval & x)
inline

Returns \(\atan([x])\).

Parameters
xinterval value
Returns
the operation result
134 {
135 Interval y = gaol::atan(x);
136 gaol::round_upward();
137 return y;
138 }

◆ atan2() [1/4]

Interval codac2::atan2 ( const Interval & y,
const Interval & x )
inline

Returns \(\mathrm{arctan2}([y],[x])\).

Parameters
yinterval value
xinterval value
Returns
the operation result
141 {
142 if(y.is_empty() || x.is_empty())
143 return Interval::empty();
144
145 // We handle the special case x=[0,0] separately
146 else if(x == Interval::zero())
147 {
148 if(y.lb() >= 0)
149 {
150 if(y.ub() == 0)
151 return Interval::empty(); // atan2(0,0) is undefined
152 else
153 return Interval::half_pi();
154 }
155
156 else if(y.ub() <= 0)
157 return -Interval::half_pi();
158
159 else
160 return Interval(-1,1)*Interval::half_pi();
161 }
162
163 else if(x.lb() >= 0)
164 return atan(y/x); // now, x.ub()>0 -> atan does not give an empty set
165
166 else if(x.ub() <= 0)
167 {
168 if(y.lb() >= 0)
169 return atan(y/x) + Interval::pi(); // x.lb()<0
170 else if(y.ub() < 0)
171 return atan(y/x) - Interval::pi();
172 else
173 return Interval(-1,1)*Interval::pi();
174 }
175
176 else
177 {
178 if(y.lb() >= 0)
179 return atan(y/x.ub()) | (atan(y/x.lb()) + Interval::pi());
180
181 else if(y.ub() <= 0)
182 {
183 if(x.lb() != -oo)
184 {
185 if(x.ub() != oo)
186 return (atan(y/x.lb())-Interval::pi()) | atan(y/x.ub());
187 else
188 return (atan(y/x.lb())-Interval::pi()) | Interval::zero();
189 }
190
191 else
192 {
193 if(x.ub() != oo)
194 return (-Interval::pi()) | atan(y/x.ub());
195 else
196 return -Interval::pi() | Interval::zero();
197 }
198 }
199
200 else
201 return Interval(-1,1)*Interval::pi();
202 }
203 }
static Interval pi()
Provides an interval for .
Definition codac2_Interval_impl.h:578
static Interval zero()
Provides an interval for .
Definition codac2_Interval_impl.h:563
static Interval half_pi()
Provides an interval for .
Definition codac2_Interval_impl.h:573
Interval atan(const Interval &x)
Returns .
Definition codac2_Interval_operations_impl.h:133

◆ cosh() [1/2]

Interval codac2::cosh ( const Interval & x)
inline

Returns \(\cosh([x])\).

Parameters
xinterval value
Returns
the operation result
206 {
207 Interval y;
208 if(x.is_unbounded())
209 y = Interval(gaol::cosh(x).left(),oo);
210 else
211 y = gaol::cosh(x);
212 gaol::round_upward();
213 return y;
214 }
bool is_unbounded() const
Tests if one of the bounds of this is infinite.
Definition codac2_Interval_impl.h:210

◆ sinh() [1/2]

Interval codac2::sinh ( const Interval & x)
inline

Returns \(\sinh([x])\).

Parameters
xinterval value
Returns
the operation result
217 {
218 Interval y = gaol::sinh(x);
219 gaol::round_upward();
220 return y;
221 }

◆ tanh() [1/2]

Interval codac2::tanh ( const Interval & x)
inline

Returns \(\tanh([x])\).

Parameters
xinterval value
Returns
the operation result
224 {
225 Interval y = gaol::tanh(x);
226 gaol::round_upward();
227 return y;
228 }

◆ acosh() [1/2]

Interval codac2::acosh ( const Interval & x)
inline

Returns \(\acosh([x])\).

Parameters
xinterval value
Returns
the operation result
231 {
232 Interval y = gaol::acosh(x);
233 gaol::round_upward();
234 return y;
235 }

◆ asinh() [1/2]

Interval codac2::asinh ( const Interval & x)
inline

Returns \(\asinh([x])\).

Parameters
xinterval value
Returns
the operation result
238 {
239 if(x.is_empty())
240 return Interval::empty();
241
242 else if(x.lb() >= 0)
243 return gaol::asinh(x);
244
245 else if(x.ub() <= 0)
246 return -gaol::asinh(-x);
247
248 else
249 return {
250 -gaol::asinh(gaol::interval(0,-x.lb())).right(),
251 gaol::asinh(gaol::interval(0,x.ub())).right()
252 };
253
254 // no round_upward?
255 }

◆ atanh() [1/2]

Interval codac2::atanh ( const Interval & x)
inline

Returns \(\atanh([x])\).

Parameters
xinterval value
Returns
the operation result
258 {
259 Interval y = gaol::atanh(x);
260 gaol::round_upward();
261 return y;
262 }

◆ abs() [1/3]

Interval codac2::abs ( const Interval & x)
inline

Returns \(\mid[x]\mid = \left\{\mid x \mid, x\in[x]\right\}\).

Parameters
xinterval value
Returns
the operation result
265 {
266 return gaol::abs(x);
267 }

◆ min() [1/4]

Interval codac2::min ( const Interval & x,
const Interval & y )
inline

Returns \(\min([x],[y])=\left\{\min(x,y), x\in[x], y\in[y]\right\}\).

Parameters
xinterval value
yinterval value
Returns
the operation result
270 {
271 return gaol::min(x,y);
272 }

◆ max() [1/4]

Interval codac2::max ( const Interval & x,
const Interval & y )
inline

Returns \(\max([x],[y])=\left\{\max(x,y), x\in[x], y\in[y]\right\}\).

Parameters
xinterval value
yinterval value
Returns
the operation result
275 {
276 return gaol::max(x,y);
277 }

◆ sign()

Interval codac2::sign ( const Interval & x)
inline

Returns \(\sign([x])=\left[\left\{\sign(x), x\in[x]\right\}\right]\).

Note
By convention, \( 0\in[x] \Longrightarrow \sign([x])=[-1,1]\).
Parameters
xinterval value
Returns
the operation result
280 {
281 return x.ub() < 0 ? -1. : x.lb() > 0 ? 1. : Interval(-1.,1.);
282 }

◆ integer()

Interval codac2::integer ( const Interval & x)
inline

Returns the largest integer interval included in \([x]\).

The result is a subset interval of x.

Parameters
xinterval value
Returns
the operation result
285 {
286 return gaol::integer(x);
287 }

◆ floor() [1/2]

Interval codac2::floor ( const Interval & x)
inline

Returns floor of \([x]\).

Parameters
xinterval value
Returns
the operation result
290 {
291 return gaol::floor(x);
292 }

◆ ceil() [1/2]

Interval codac2::ceil ( const Interval & x)
inline

Returns ceil of \([x]\).

Parameters
xinterval value
Returns
the operation result
295 {
296 return gaol::ceil(x);
297 }

◆ chi()

Interval codac2::chi ( const Interval & x,
const Interval & y,
const Interval & z )
inline

Return \([y]\) if \(x^+ <= 0\), \([z]\) if \(x^- > 0\), \([y]\sqcup[z]\) else.

Parameters
xinterval value
yinterval value
zinterval value
Returns
the operation result
300 {
301 if(x.ub() <= 0)
302 return y;
303
304 else if(x.lb() > 0)
305 return z;
306
307 else
308 return y | z;
309 }

◆ operator<<() [4/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const IntervalMatrix & x )
inline

Stream output operator for IntervalMatrix objects.

Parameters
osThe output stream to write to.
xThe interval matrix whose contents are to be printed.
Returns
A reference to the modified output stream.
35 {
36 if(x.is_empty())
37 return os << "[ empty " << x.rows() << "x" << x.cols() << " mat ]";
38
39 else
40 {
41 os << x.format(codac_matrix_fmt());
42 return os;
43 }
44 }
Eigen::IOFormat codac_matrix_fmt()
Provides an Eigen IOFormat for formatting matrices.
Definition codac2_matrices.h:232

◆ operator<<() [5/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const IntervalRow & x )
inline

Stream output operator for IntervalRow objects.

Parameters
osThe output stream to write to.
xThe interval row whose contents are to be printed.
Returns
A reference to the modified output stream.
34 {
35 if(x.is_empty())
36 return os << "[ empty row ]";
37
38 else
39 {
40 os << x.format(codac_row_fmt());
41 return os;
42 }
43 }
Eigen::IOFormat codac_row_fmt()
Provides an Eigen IOFormat for formatting row vectors.
Definition codac2_matrices.h:205

◆ operator<<() [6/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const IntervalVector & x )
inline

Stream output operator for IntervalVector objects.

Parameters
osThe output stream to write to.
xThe interval vector whose contents are to be printed.
Returns
A reference to the modified output stream.
35 {
36 if(x.is_empty())
37 return os << "[ empty vector ]";
38
39 else
40 {
41 os << x.format(codac_vector_fmt());
42 return os;
43 }
44 }
Eigen::IOFormat codac_vector_fmt()
Provides an Eigen IOFormat for formatting column vectors.
Definition codac2_matrices.h:218

◆ operator&() [2/3]

ConvexPolygon codac2::operator& ( const ConvexPolygon & p1,
const ConvexPolygon & p2 )

Computes the intersection of two convex polygons.

The result is a new ConvexPolygon enclosing the intersection area. If the polygons do not intersect, the result is an empty convex polygon.

Parameters
p1The first convex polygon.
p2The second convex polygon.
Returns
A ConvexPolygon enclosing the intersection area.

◆ operator<<() [7/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const OrientationInterval & x )
inline

Streams out a OrientationInterval.

Parameters
osThe stream to be updated
xThe orientation interval to stream out
Returns
A Reference to the updated stream
45 {
46 if(x == OrientationInterval::EMPTY)
47 os << "[ empty ]";
48 else if(x == OrientationInterval::COLINEAR)
49 os << "[ col ]";
50 else if(x == OrientationInterval::CLOCKWISE)
51 os << "[ cw ]";
52 else if(x == OrientationInterval::COUNTERCLOCKWISE)
53 os << "[ ccw ]";
54 else if(x == (OrientationInterval::COLINEAR | OrientationInterval::CLOCKWISE))
55 os << "[ col, cw ]";
56 else if(x == (OrientationInterval::COLINEAR | OrientationInterval::COUNTERCLOCKWISE))
57 os << "[ col, ccw ]";
59 os << "[ col, cw, ccw ]";
60 return os;
61 }
@ UNKNOWN
Definition codac2_geometry.h:31

◆ orientation()

OrientationInterval codac2::orientation ( const IntervalVector & p1,
const IntervalVector & p2,
const IntervalVector & p3 )

Computes the orientation of an ordered triplet of 2D points.

Determines whether the oriented angle \(\widehat{p_1 p_2 p_3}\) is positive (counterclockwise), negative (clockwise), or if the points are colinear (flat or 0 angle). Depending on floating point uncertainties, the test may not be able to conclude (a UNKNOWN value would then be returned).

Parameters
p1First point (2d IntervalVector) of the triplet.
p2Second point (2d IntervalVector) of the triplet (vertex of the angle).
p3Third point (2d IntervalVector) of the triplet.
Returns
An orientation of type OrientationInterval

◆ aligned()

BoolInterval codac2::aligned ( const IntervalVector & p1,
const IntervalVector & p2,
const IntervalVector & p3 )

Checks whether three 2D points are aligned (colinear).

Determines if the points lie on the same straight line using an orientation test (cross product). Depending on floating point uncertainties, the test may not be able to conclude (a UNKNOWN value would then be returned).

Parameters
p1First point (2d IntervalVector) of the triplet.
p2Second point (2d IntervalVector) of the triplet.
p3Third point (2d IntervalVector) of the triplet.
Returns
A BooleanInterval

◆ convex_hull()

std::vector< IntervalVector > codac2::convex_hull ( std::vector< IntervalVector > pts)

Computes the convex hull of a set of 2d points.

Given a set of 2d points enclosed in tiny boxes, the function computes their convex hull. The method is based on a Graham scan algorithm. The output list of the algorithm is a subset of the input list, with same uncertainties and a possible different order.

Parameters
pts2d points in any order
Returns
Points on the convex hull in counterclockwise order.

◆ rotate() [1/2]

IntervalVector codac2::rotate ( const IntervalVector & p,
const Interval & a,
const IntervalVector & c = Vector::zero(2) )

Rotates a 2D interval vector around a point by a given angle.

This function performs a 2D rotation of the input interval vector p by an interval angle a, optionally around a center of rotation c (default is the origin). This operation propagates uncertainties through interval arithmetic.

Parameters
pThe input 2D interval vector to be rotated.
aThe rotation angle as an interval (in radians).
cThe center of rotation as an interval vector (default is the origin (0,0)).
Returns
The rotated interval vector.

◆ rotate() [2/2]

template<typename T>
T codac2::rotate ( const T & x,
const Interval & a,
const IntervalVector & c = Vector::zero(2) )
inline

Rotates a 2D object (Segment, Polygon, etc) around a point by a given angle.

This function performs a 2D rotation of the input object x by an interval angle a, optionally around a center of rotation c (default is the origin). This operation propagates uncertainties through interval arithmetic.

Parameters
xThe input 2D object to be rotated.
aThe rotation angle as an interval (in radians).
cThe center of rotation as an interval vector (default is the origin (0,0)).
Returns
The rotated object.
136 {
137 T x_(x);
138 for(auto& ei : x_)
139 ei = rotate(ei,a,c);
140 return x_;
141 }
IntervalVector rotate(const IntervalVector &p, const Interval &a, const IntervalVector &c=Vector::zero(2))
Rotates a 2D interval vector around a point by a given angle.

◆ operator<<() [8/12]

std::ostream & codac2::operator<< ( std::ostream & str,
const Polygon & p )

Stream output operator for Polygon.

Parameters
strOutput stream.
pThe polygon to print.
Returns
The output stream with the polygon information.

◆ operator&() [3/3]

IntervalVector codac2::operator& ( const Segment & e1,
const Segment & e2 )

Computes the intersection of two segments.

If the segments do not intersect, an empty IntervalVector is returned. If the segments are colinear, the set of intersection points is returned as a box.

Parameters
e1The first segment.
e2The second segment.
Returns
An IntervalVector enclosing the intersection point.

◆ proj_intersection()

IntervalVector codac2::proj_intersection ( const Segment & e1,
const Segment & e2 )

Computes the projected intersection of two segments.

This corresponds to the intersection of the two lines related to the two segments. Therefore, the intersection point may not belong to the segments.

If the segments are parallel but not colinear, an empty IntervalVector is returned. If the segments are colinear, the set of intersection points is returned as a box.

Parameters
e1The first segment.
e2The second segment.
Returns
An IntervalVector enclosing the intersection point.

◆ colinear()

BoolInterval codac2::colinear ( const Segment & e1,
const Segment & e2 )

Checks if two segments are colinear.

Parameters
e1The first segment.
e2The second segment.
Returns
A BoolInterval indicating possible colinearity.

◆ operator<<() [9/12]

std::ostream & codac2::operator<< ( std::ostream & str,
const Segment & e )

Stream output operator for Segment.

Parameters
strOutput stream.
eThe segment to print.
Returns
The output stream with the segment information.

◆ inverse_correction()

template<LeftOrRightInv O = LEFT_INV, typename OtherDerived, typename OtherDerived_>
IntervalMatrix codac2::inverse_correction ( const Eigen::MatrixBase< OtherDerived > & A,
const Eigen::MatrixBase< OtherDerived_ > & B )
inline

Correct the approximation of the inverse \(\mathbf{B}\approx\mathbf{A}^{-1}\) of a square matrix \(\mathbf{A}\) by providing a reliable enclosure \([\mathbf{A}^{-1}]\).

Precondition
\(\mathbf{A}\) and \(\mathbf{B}\) are square matrices, possibly interval matrices.
Template Parameters
OIf LEFT_INV, use the inverse of \(\mathbf{BA}\) (otherwise use the inverse of \(\mathbf{AB}\), left inverse is normally better). In Python/Matlab, this template parameter is provided as a last boolean argument, and is left_inv = True by default.
Parameters
AA matrix expression, possibly interval.
BAn (almost punctual) approximation of its inverse.
Returns
The enclosure of the inverse.
36 {
37 assert_release(A.is_squared());
38 assert_release(B.is_squared());
39
40 auto A_ = A.template cast<Interval>();
41 auto B_ = B.template cast<Interval>();
42
43 Index N = A_.rows();
44 assert_release(N==B_.rows());
45
46 auto Id = IntervalMatrix::Identity(N,N);
47 auto erMat = [&]() { if constexpr(O == LEFT_INV) return -B_*A_+Id; else return -A_*B_+Id; }();
48
49 double mrad=0.0;
51 IntervalMatrix Ep = Id+erMat*(Id+E);
52 /* one could use several iterations here, either
53 using mrad, or directly */
54
55 auto res = (O == LEFT_INV) ? IntervalMatrix(Ep*B_) : IntervalMatrix(B_*Ep);
56
57 // We want the presence of non-invertible matrices to
58 // result in coefficients of the form [-oo,+oo].
59 if (mrad==oo) {
60 for (Index c=0;c<N;c++) {
61 for (Index r=0;r<N;r++) {
62 if (Ep(r,c).is_unbounded()) {
63 for (Index k=0;k<N;k++) {
64 if constexpr(O == LEFT_INV)
65 res(r,k) = Interval();
66 else
67 res(k,c) = Interval();
68 }
69 }
70 }
71 }
72 }
73 return res;
74 }
bool is_unbounded() const
Checks if the interval matrix contains any unbounded intervals.
Definition codac2_MatrixBase_addons_IntervalMatrixBase.h:283
IntervalMatrix infinite_sum_enclosure(const IntervalMatrix &A, double &mrad)
Compute an upper bound of , with a matrix of intervals as an "error term" (uses only bounds on coeff...
Eigen::Matrix< Interval,-1,-1 > IntervalMatrix
Alias for a dynamic-size matrix of intervals.
Definition codac2_IntervalMatrix.h:25

◆ inverse_enclosure()

template<typename OtherDerived>
IntervalMatrix codac2::inverse_enclosure ( const Eigen::MatrixBase< OtherDerived > & A)
inline

Enclosure of the inverse of a (non-singular) matrix expression, possibly an interval matrix.

Precondition
\(\mathbf{A}\) is a square matrix.
Parameters
AA matrix expression, possibly interval.
Returns
The enclosure of the inverse. Can have \((-\infty,\infty)\) coefficients if \(\mathbf{A}\) is singular or almost singular, if the inversion "failed".
88 {
89 assert_release(A.is_squared());
90 Index N = A.rows();
91
92 if constexpr(std::is_same_v<typename OtherDerived::Scalar,Interval>)
94 (A.mid()).fullPivLu().solve(Matrix::Identity(N,N)));
95
96 else
98 A.fullPivLu().solve(Matrix::Identity(N,N)));
99 }
IntervalMatrix inverse_correction(const Eigen::MatrixBase< OtherDerived > &A, const Eigen::MatrixBase< OtherDerived_ > &B)
Correct the approximation of the inverse of a square matrix by providing a reliable enclosure .
Definition codac2_inversion.h:35

◆ infinite_sum_enclosure()

IntervalMatrix codac2::infinite_sum_enclosure ( const IntervalMatrix & A,
double & mrad )

Compute an upper bound of \(\left([\mathbf{A}]+[\mathbf{A}]^2+[\mathbf{A}]^3+\dots\right)\), with \([\mathbf{A}]\) a matrix of intervals as an "error term" (uses only bounds on coefficients).

The function also returns mrad, which gives an idea of the magnification of the matrix during calculation. In particular, if mrad = \(\infty\), then the inversion calculation (e.g., performed by Eigen) has somehow failed and some coefficients of the output interval matrix are \([-\infty,\infty]\).

Precondition
\([\mathbf{A}]\) is a square matrix.
Parameters
AA matrix of intervals (supposed around \(\mathbf{0}\)).
mradThe maximum radius of the result added (output argument).
Returns
The sum enclosure. May be unbounded.

◆ abs() [2/3]

template<typename OtherDerived>
auto codac2::abs ( const Eigen::MatrixBase< OtherDerived > & x)
inline

Compute the element-wise absolute value of a matrix.

This function takes an Eigen matrix expression and returns a matrix where each element is replaced by its absolute value.

For scalar type double, the standard library fabs() is used. For other scalar types, it uses the generic abs() function.

Parameters
xInput matrix expression.
Returns
A new Eigen matrix with the absolute values of the elements of x.
127 {
128 using M = Eigen::MatrixBase<OtherDerived>;
129 Eigen::Matrix<typename M::Scalar,M::RowsAtCompileTime,M::ColsAtCompileTime> a(x.rows(),x.cols());
130
131 for(Index i = 0 ; i < x.rows() ; i++)
132 for(Index j = 0 ; j < x.cols() ; j++)
133 {
134 if constexpr(std::is_same_v<typename M::Scalar,double>)
135 a(i,j) = fabs(x(i,j));
136 else
137 a(i,j) = abs(x(i,j));
138 }
139
140 return a;
141 }
Interval abs(const Interval &x)
Returns .
Definition codac2_Interval_operations_impl.h:264

◆ floor() [2/2]

template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto codac2::floor ( const Eigen::MatrixBase< OtherDerived > & x)
inline

Compute the element-wise floor of a matrix.

This function returns a matrix where each element is replaced by the largest integer not greater than that element.

Disabled for interval matrices.

Parameters
xInput matrix expression.
Returns
A new Eigen matrix with floored elements.
157 {
158 return x.array().floor().matrix();
159 }

◆ ceil() [2/2]

template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto codac2::ceil ( const Eigen::MatrixBase< OtherDerived > & x)
inline

Compute the element-wise ceiling of a matrix.

This function returns a matrix where each element is replaced by the smallest integer not less than that element.

Disabled for interval matrices.

Parameters
xInput matrix expression.
Returns
A new Eigen matrix with ceiled elements.
175 {
176 return x.array().ceil().matrix();
177 }

◆ round()

template<typename OtherDerived>
requires (!Eigen::IsIntervalDomain<typename OtherDerived::Scalar>)
auto codac2::round ( const Eigen::MatrixBase< OtherDerived > & x)
inline

Compute the element-wise rounding of a matrix.

This function returns a matrix where each element is replaced by the nearest integer to that element.

Disabled for interval matrices.

Parameters
xInput matrix expression.
Returns
A new Eigen matrix with rounded elements.
193 {
194 return x.array().round().matrix();
195 }

◆ codac_row_fmt()

Eigen::IOFormat codac2::codac_row_fmt ( )
inline

Provides an Eigen IOFormat for formatting row vectors.

This format prints elements separated by spaces, with brackets around the entire row vector.

Returns
An Eigen::IOFormat configured for row vector formatting.
206 {
207 return Eigen::IOFormat(Eigen::StreamPrecision, Eigen::DontAlignCols, " ", "", "", "", "[ ", " ]");
208 }

◆ codac_vector_fmt()

Eigen::IOFormat codac2::codac_vector_fmt ( )
inline

Provides an Eigen IOFormat for formatting column vectors.

This format prints elements separated by semicolons, with brackets around the entire vector.

Returns
An Eigen::IOFormat configured for column vector formatting.
219 {
220 return Eigen::IOFormat(Eigen::StreamPrecision, Eigen::DontAlignCols, "", " ; ", "", "", "[ ", " ]");
221 }

◆ codac_matrix_fmt()

Eigen::IOFormat codac2::codac_matrix_fmt ( )
inline

Provides an Eigen IOFormat for formatting matrices.

This format prints elements separated by commas, rows separated by new lines, and brackets around the entire matrix.

Returns
An Eigen::IOFormat configured for matrix formatting.
233 {
234 return Eigen::IOFormat(Eigen::StreamPrecision, 0, " , ", "\n", "[ ", " ]", "[", "]");
235 }

◆ operator<<() [10/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const Matrix & x )
inline

Stream output operator for Matrix objects.

Parameters
osThe output stream to write to.
xThe matrix whose contents are to be printed.
Returns
A reference to the modified output stream.
36 {
37 os << x.format(codac_matrix_fmt());
38 return os;
39 }

◆ operator<<() [11/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const Row & x )
inline

Stream output operator for Row objects.

Parameters
osThe output stream to write to.
xThe row whose contents are to be printed.
Returns
A reference to the modified output stream.
34 {
35 os << x.format(codac_row_fmt());
36 return os;
37 }

◆ operator<<() [12/12]

std::ostream & codac2::operator<< ( std::ostream & os,
const Vector & x )
inline

Stream output operator for Vector objects.

Parameters
osThe output stream to write to.
xThe vector whose contents are to be printed.
Returns
A reference to the modified output stream.
34 {
35 os << x.format(codac_vector_fmt());
36 return os;
37 }

◆ serialize() [1/4]

template<typename T>
requires (std::is_trivially_copyable_v<T>)
void codac2::serialize ( std::ostream & f,
const T & x )
inline

Writes the binary representation of a trivially copyable object to the given output stream.

Binary structure:
[raw memory of x]

Parameters
foutput stream
xobject/variable to be serialized
33 {
34 f.write(reinterpret_cast<const char*>(&x), sizeof(T));
35 }

◆ deserialize() [1/4]

template<typename T>
requires (std::is_trivially_copyable_v<T>)
void codac2::deserialize ( std::istream & f,
T & x )
inline

Reads the binary representation of a trivially copyable object from the given input stream.

Binary structure:
[raw memory of x]

Parameters
finput stream
xobject to be deserialized
50 {
51 f.read(reinterpret_cast<char*>(&x), sizeof(T));
52 }

◆ serialize() [2/4]

void codac2::serialize ( std::ostream & f,
const Interval & x )
inline

Writes the binary representation of an Interval object to the given output stream.

Interval binary structure:
[double_lb][double_ub]

Parameters
foutput stream
xInterval object to be serialized
67 {
68 serialize(f,x.lb());
69 serialize(f,x.ub());
70 }
void serialize(std::ostream &f, const T &x)
Writes the binary representation of a trivially copyable object to the given output stream.
Definition codac2_serialization.h:32

◆ deserialize() [2/4]

void codac2::deserialize ( std::istream & f,
Interval & x )
inline

Creates an Interval object from the binary representation given in an input stream.

Interval binary structure:
[double_lb][double_ub]

Parameters
finput stream
xInterval object to be deserialized
83 {
84 double lb, ub;
85 deserialize(f,lb);
86 deserialize(f,ub);
87 x = Interval(lb,ub);
88 }
auto lb() const
Returns a matrix containing the lower bounds of each interval element.
Definition codac2_MatrixBase_addons_IntervalMatrixBase.h:91
auto ub() const
Returns a matrix containing the upper bounds of each interval element.
Definition codac2_MatrixBase_addons_IntervalMatrixBase.h:103
void deserialize(std::istream &f, T &x)
Reads the binary representation of a trivially copyable object from the given input stream.
Definition codac2_serialization.h:49

◆ serialize() [3/4]

template<typename T, int R = -1, int C = -1>
void codac2::serialize ( std::ostream & f,
const Eigen::Matrix< T, R, C > & x )
inline

Writes the binary representation of an Eigen::Matrix to the given output stream. The structure can be a matrix, a row, a vector, with double or Interval components.

Matrix binary structure:
[Index_rows][Index_cols][element_00][element_01]...[element_rc]

Parameters
foutput stream
xmatrix structure to be serialized
105 {
106 Index r = x.rows(), c = x.cols();
107 f.write(reinterpret_cast<const char*>(&r), sizeof(Index));
108 f.write(reinterpret_cast<const char*>(&c), sizeof(Index));
109 for(Index i = 0 ; i < r ; i++)
110 for(Index j = 0 ; j < c ; j++)
111 serialize(f,x(i,j));
112 }

◆ deserialize() [3/4]

template<typename T, int R = -1, int C = -1>
void codac2::deserialize ( std::istream & f,
Eigen::Matrix< T, R, C > & x )
inline

Reads the binary representation of an Eigen::Matrix from the given input stream. The structure can be a matrix, a row, a vector, with double or Interval components.

Matrix binary structure:
[Index_rows][Index_cols][element_00][element_01]...[element_rc]

Parameters
finput stream
xmatrix structure to be deserialized
127 {
128 Index r, c;
129 f.read(reinterpret_cast<char*>(&r), sizeof(Index));
130 f.read(reinterpret_cast<char*>(&c), sizeof(Index));
131
132 if constexpr(R == -1 && C == -1)
133 x.resize(r,c);
134 else if constexpr(R == -1 || C == -1)
135 x.resize(std::max(r,c));
136 else
137 assert_release(R == r && C == c);
138
139 for(Index i = 0 ; i < r ; i++)
140 for(Index j = 0 ; j < c ; j++)
141 deserialize(f,x(i,j));
142 }

◆ serialize() [4/4]

template<typename T>
void codac2::serialize ( std::ostream & f,
const SampledTraj< T > & x )
inline

Writes the binary representation of a SampledTraj object to the given output stream.

SampledTraj binary structure:
[nb_samples][t0][x0][t1][x1]...[tn][xn]

Parameters
foutput stream
xSampledTraj object to be serialized
158 {
159 Index size = x.nb_samples();
160 f.write(reinterpret_cast<const char*>(&size), sizeof(size));
161
162 for(const auto& [ti, xi] : x)
163 {
164 serialize(f, ti);
165 serialize(f, xi);
166 }
167 }

◆ deserialize() [4/4]

template<typename T>
void codac2::deserialize ( std::istream & f,
SampledTraj< T > & x )
inline

Reads the binary representation of a SampledTraj object from the given input stream.

SampledTraj binary structure:
[nb_samples][t0][x0][t1][x1]...[tn][xn]

Parameters
finput stream
xSampledTraj object to be deserialized
181 {
182 Index size;
183 x.clear();
184 f.read(reinterpret_cast<char*>(&size), sizeof(size));
185
186 for(Index i = 0 ; i < size ; i++)
187 {
188 double ti; deserialize(f,ti);
189 T xi; deserialize(f,xi);
190 x.set(xi,ti);
191 }
192 }

◆ operator+() [6/9]

template<typename T>
const SampledTraj< T > & codac2::operator+ ( const SampledTraj< T > & x1)
inline

\(x1(\cdot)\)

Parameters
x1
Returns
trajectory output
107 {
108 return x1;
109 }

◆ operator+() [7/9]

template<typename T>
SampledTraj< T > codac2::operator+ ( const SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

\(x_1(\cdot)+x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator+() [8/9]

template<typename T, typename Q>
SampledTraj< T > codac2::operator+ ( const SampledTraj< T > & x1,
const Q & x2 )
inline

\(x_1(\cdot)+x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator+() [9/9]

template<typename T, typename Q>
SampledTraj< T > codac2::operator+ ( const Q & x1,
const SampledTraj< T > & x2 )
inline

\(x+x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator+=() [1/2]

template<typename T, typename Q>
SampledTraj< T > & codac2::operator+= ( SampledTraj< T > & x1,
const Q & x2 )
inline

Operates +=.

Parameters
x1
x2
Returns
updated output

◆ operator+=() [2/2]

template<typename T>
SampledTraj< T > & codac2::operator+= ( SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

Operates +=.

Parameters
x1
x2
Returns
updated output

◆ operator-() [4/7]

template<typename T>
SampledTraj< T > codac2::operator- ( const SampledTraj< T > & x1)
inline

\(-x_1(\cdot)\)

Parameters
x1
Returns
trajectory output
163 {
164 return -1.*x1;
165 }

◆ operator-() [5/7]

template<typename T>
SampledTraj< T > codac2::operator- ( const SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

\(x_1(\cdot)-x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator-() [6/7]

template<typename T, typename Q>
SampledTraj< T > codac2::operator- ( const SampledTraj< T > & x1,
const Q & x2 )
inline

\(x_1(\cdot)-x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator-() [7/7]

template<typename T, typename Q>
SampledTraj< T > codac2::operator- ( const Q & x1,
const SampledTraj< T > & x2 )
inline

\(x-x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator-=() [1/2]

template<typename T, typename Q>
SampledTraj< T > & codac2::operator-= ( SampledTraj< T > & x1,
const Q & x2 )
inline

Operates -=.

Parameters
x1
x2
Returns
updated output

◆ operator-=() [2/2]

template<typename T>
SampledTraj< T > & codac2::operator-= ( SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

Operates -=.

Parameters
x1
x2
Returns
updated output

◆ operator*() [4/9]

template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > codac2::operator* ( double x1,
const SampledTraj< T > & x2 )
inline

\(x_1\cdot x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*() [5/9]

template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > codac2::operator* ( const SampledTraj< T > & x1,
double x2 )
inline

\(x_1(\cdot)\cdot x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*() [6/9]

template<typename T>
SampledTraj< T > codac2::operator* ( const SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

\(x_1(\cdot)\cdot x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*() [7/9]

template<typename T, typename Q>
SampledTraj< T > codac2::operator* ( const SampledTraj< T > & x1,
const Q & x2 )
inline

\(x_1(\cdot)\cdot x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*() [8/9]

template<typename T, typename Q>
SampledTraj< T > codac2::operator* ( const Q & x1,
const SampledTraj< T > & x2 )
inline

\(x\cdot x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*() [9/9]

SampledTraj< Vector > codac2::operator* ( const SampledTraj< Matrix > & x1,
const SampledTraj< Vector > & x2 )
inline

\(x_1(\cdot)\cdot x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator*=() [1/2]

template<typename T, typename Q>
SampledTraj< T > & codac2::operator*= ( SampledTraj< T > & x1,
const Q & x2 )
inline

Operates *=.

Parameters
x1
x2
Returns
updated output

◆ operator*=() [2/2]

template<typename T>
SampledTraj< T > & codac2::operator*= ( SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

Operates *=.

Parameters
x1
x2
Returns
updated output

◆ operator/() [4/7]

template<typename T>
requires (!std::is_same_v<T,double>)
SampledTraj< T > codac2::operator/ ( const SampledTraj< T > & x1,
double x2 )
inline

\(x_2(\cdot)/x_1\)

Parameters
x1
x2
Returns
trajectory output

◆ operator/() [5/7]

template<typename T>
SampledTraj< T > codac2::operator/ ( const SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

\(x_1(\cdot)/x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator/() [6/7]

template<typename T, typename Q>
SampledTraj< T > codac2::operator/ ( const SampledTraj< T > & x1,
const Q & x2 )
inline

\(x_1(\cdot)/x_2\)

Parameters
x1
x2
Returns
trajectory output

◆ operator/() [7/7]

template<typename T, typename Q>
SampledTraj< T > codac2::operator/ ( const Q & x1,
const SampledTraj< T > & x2 )
inline

\(x/x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ operator/=() [1/2]

template<typename T, typename Q>
SampledTraj< T > & codac2::operator/= ( SampledTraj< T > & x1,
const Q & x2 )
inline

Operates /=.

Parameters
x1
x2
Returns
updated output

◆ operator/=() [2/2]

template<typename T>
SampledTraj< T > & codac2::operator/= ( SampledTraj< T > & x1,
const SampledTraj< T > & x2 )
inline

Operates /=.

Parameters
x1
x2
Returns
updated output

◆ sqr() [2/2]

SampledTraj< double > codac2::sqr ( const SampledTraj< double > & x1)

\(x^2(\cdot)\)

Parameters
x1
Returns
trajectory output

◆ sqrt() [2/2]

SampledTraj< double > codac2::sqrt ( const SampledTraj< double > & x1)

\(\sqrt{x_1(\cdot)}\)

Parameters
x1
Returns
trajectory output

◆ pow() [4/5]

SampledTraj< double > codac2::pow ( const SampledTraj< double > & x1,
int x2 )

\(x^x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ pow() [5/5]

SampledTraj< double > codac2::pow ( const SampledTraj< double > & x1,
double x2 )

\(x^x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ root() [2/2]

SampledTraj< double > codac2::root ( const SampledTraj< double > & x1,
int x2 )

\(x^x_2(\cdot)\)

Parameters
x1
x2
Returns
trajectory output

◆ exp() [2/2]

SampledTraj< double > codac2::exp ( const SampledTraj< double > & x1)

\(\exp(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ log() [2/2]

SampledTraj< double > codac2::log ( const SampledTraj< double > & x1)

\(\log(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ cos() [2/2]

SampledTraj< double > codac2::cos ( const SampledTraj< double > & x1)

\(\cos(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ sin() [2/2]

SampledTraj< double > codac2::sin ( const SampledTraj< double > & x1)

\(\sin(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ tan() [2/2]

SampledTraj< double > codac2::tan ( const SampledTraj< double > & x1)

\(\tan(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ acos() [2/2]

SampledTraj< double > codac2::acos ( const SampledTraj< double > & x1)

\(\arccos(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ asin() [2/2]

SampledTraj< double > codac2::asin ( const SampledTraj< double > & x1)

\(\arcsin(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ atan() [2/2]

SampledTraj< double > codac2::atan ( const SampledTraj< double > & x1)

\(\arctan(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ atan2() [2/4]

SampledTraj< double > codac2::atan2 ( const SampledTraj< double > & x1,
const SampledTraj< double > & x2 )

\(\mathrm{arctan2}(x_1(\cdot),x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ atan2() [3/4]

SampledTraj< double > codac2::atan2 ( const SampledTraj< double > & x1,
double x2 )

\(\mathrm{arctan2}(x_1(\cdot),x_2)\)

Parameters
x1
x2
Returns
trajectory output

◆ atan2() [4/4]

SampledTraj< double > codac2::atan2 ( double x1,
const SampledTraj< double > & x2 )

\(\mathrm{arctan2}(x_1, x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ cosh() [2/2]

SampledTraj< double > codac2::cosh ( const SampledTraj< double > & x1)

\(\cosh(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ sinh() [2/2]

SampledTraj< double > codac2::sinh ( const SampledTraj< double > & x1)

\(\sinh(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ tanh() [2/2]

SampledTraj< double > codac2::tanh ( const SampledTraj< double > & x1)

\(\tanh(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ acosh() [2/2]

SampledTraj< double > codac2::acosh ( const SampledTraj< double > & x1)

\(\mathrm{arccosh}(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ asinh() [2/2]

SampledTraj< double > codac2::asinh ( const SampledTraj< double > & x1)

\(\mathrm{arcsinh}(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ atanh() [2/2]

SampledTraj< double > codac2::atanh ( const SampledTraj< double > & x1)

\(\mathrm{arctanh}(x_1(\cdot))\)

Parameters
x1
Returns
trajectory output

◆ abs() [3/3]

SampledTraj< double > codac2::abs ( const SampledTraj< double > & x1)

\(\mid x_1(\cdot)\mid\)

Parameters
x1
Returns
trajectory output

◆ min() [2/4]

SampledTraj< double > codac2::min ( const SampledTraj< double > & x1,
const SampledTraj< double > & x2 )

\(\min(x_1(\cdot),x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ min() [3/4]

SampledTraj< double > codac2::min ( const SampledTraj< double > & x1,
double x2 )

\(\min(x_1(\cdot),x_2)\)

Parameters
x1
x2
Returns
trajectory output

◆ min() [4/4]

SampledTraj< double > codac2::min ( double x1,
const SampledTraj< double > & x2 )

\(\min(x_1, x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ max() [2/4]

SampledTraj< double > codac2::max ( const SampledTraj< double > & x1,
const SampledTraj< double > & x2 )

\(\max(x_1(\cdot),x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ max() [3/4]

SampledTraj< double > codac2::max ( const SampledTraj< double > & x1,
double x2 )

\(\max(x_1(\cdot),x_2)\)

Parameters
x1
x2
Returns
trajectory output

◆ max() [4/4]

SampledTraj< double > codac2::max ( double x1,
const SampledTraj< double > & x2 )

\(\max(x_1, x_2(\cdot))\)

Parameters
x1
x2
Returns
trajectory output

◆ to_capd() [1/3]

capd::Interval codac2::to_capd ( const codac2::Interval & x)

Converts a Codac Interval object into an CAPD Interval object.

Parameters
xthe Codac Interval to be converted
Returns
the converted CAPD Interval

◆ to_codac() [1/3]

codac2::Interval codac2::to_codac ( const capd::Interval & x)

Converts an CAPD Interval object into a Codac Interval object.

Parameters
xthe CAPD Interval to be converted
Returns
the converted Codac Interval

◆ to_capd() [2/3]

capd::IVector codac2::to_capd ( const codac2::IntervalVector & x)

Converts a Codac IntervalVector object into an CAPD IVector object.

Parameters
xthe Codac IntervalVector to be converted
Returns
the converted CAPD IVector

◆ to_codac() [2/3]

codac2::IntervalVector codac2::to_codac ( const capd::IVector & x)

Converts an CAPD IVector object into a Codac IntervalVector object.

Parameters
xthe CAPD IVector to be converted
Returns
the converted Codac IntervalVector

◆ to_capd() [3/3]

capd::IMatrix codac2::to_capd ( const codac2::IntervalMatrix & x)

Converts a Codac IntervalMatrix object into an CAPD IMatrix object.

Parameters
xthe Codac IntervalMatrix to be converted
Returns
the converted CAPD IMatrix

◆ to_codac() [3/3]

codac2::IntervalMatrix codac2::to_codac ( const capd::IMatrix & x)

Converts an CAPD IMatrix object into a Codac IntervalMatrix object.

Parameters
xthe CAPD IMatrix to be converted
Returns
the converted Codac IntervalMatrix

◆ draw_while_paving() [1/4]

void codac2::draw_while_paving ( const IntervalVector & x0,
std::shared_ptr< const CtcBase< IntervalVector > > c,
double eps,
std::shared_ptr< Figure2D > fig = nullptr )

Draws the paving while it is being computed.

Parameters
x0Initial box
cContractor to use
epsAccuracy of the paving algorithm (the undefined boxes will have their max_diam <= eps)
figFigure to draw on (optionnal, if not provided, the default figure will be used)

◆ draw_while_paving() [2/4]

void codac2::draw_while_paving ( const IntervalVector & x0,
const CtcBase< IntervalVector > & c,
double eps,
std::shared_ptr< Figure2D > fig = nullptr )

Draws the paving while it is being computed.

Parameters
x0Initial box
cContractor to use
epsAccuracy of the paving algorithm (the undefined boxes will have their max_diam <= eps)
figFigure to draw on (optionnal, if not provided, the default figure will be used)

◆ draw_while_paving() [3/4]

void codac2::draw_while_paving ( const IntervalVector & x0,
std::shared_ptr< const SepBase > s,
double eps,
std::shared_ptr< Figure2D > fig = nullptr )

Draws the paving while it is being computed.

Parameters
x0Initial box
sSeparator to use
epsAccuracy of the paving algorithm (the undefined boxes will have their max_diam <= eps)
figFigure to draw on (optionnal, if not provided, the default figure will be used)

◆ draw_while_paving() [4/4]

void codac2::draw_while_paving ( const IntervalVector & x0,
const SepBase & s,
double eps,
std::shared_ptr< Figure2D > fig = nullptr )

Draws the paving while it is being computed.

Parameters
x0Initial box
sSeparator to use
epsAccuracy of the paving algorithm (the undefined boxes will have their max_diam <= eps)
figFigure to draw on (optionnal, if not provided, the default figure will be used)