Lattice Tester Online Documentation unknown
Software Package For Testing The Uniformity Of Integral Lattices In The Real Space
Loading...
Searching...
No Matches
LatticeTester Namespace Reference

LatticeTester namespace. More...

Namespaces

namespace  CoordinateSets
 The classes FromRange, SubSets, and AddCoordinate are defined here.
 
namespace  Random
 

Classes

class  Chrono
 This class provides Chrono objects that act as stopwatches that use the system clock to compute the CPU time used by parts of a program. More...
 
class  Coordinates
 An object type that contains a set of coordinate indices, used to specify a projection. More...
 
class  FigureOfMeritDualM
 This class offers tools to calculate the same figure of merit (FOM) as FigureOfMerit, but for the m-duals of the projections. More...
 
class  FigureOfMeritM
 This class provides tools to calculate the figure of merit (FOM) \( M_{t_1,\dots,t_d}\) for any given IntLatticeExt object. More...
 
class  IntLattice
 An IntLattice object is an integral lattice, with its basis or its m-dual basis, or both. More...
 
class  IntLatticeExt
 This abstract class extends IntLattice and is a skeleton for the specialized subclasses that define specific types of lattices. More...
 
class  MRGLattice
 This subclass of IntLatticeExt defines an MRG lattice. More...
 
class  NormaBestLat
 This Normalizer class implements approximate upper bounds on the length of the shortest nonzero vector in a lattice. More...
 
class  NormaBestUpBound
 In this normalizer, the Hermite constants \(\gamma_s\) are approximated using the best upper bounds that are available. More...
 
class  NormaLaminated
 This Normalizer class implements approximate upper bounds on the length of the shortest nonzero vector in a lattice. More...
 
class  Normalizer
 This is a base class for implementing normalization constants used in figures of merit, to normalize the length of the shortest nonzero vector in either the primal or dual lattice. More...
 
class  NormaMinkHlaw
 This class implements lower bounds on the Hermite constants based on the Minkowski-Hlawka theorem [mHLA43a]. More...
 
class  NormaPalpha
 This class implements theoretical bounds on the values of \(P_{\alpha}\) for a lattice (see class Palpha). More...
 
class  NormaRogers
 This class implements upper bounds on the length of the shortest nonzero vector in a lattice, in which the Hermite constants \(\gamma_s\) are approximated by their Rogers's bounds. More...
 
class  Rank1Lattice
 This subclass of IntLatticeExt defines a general rank 1 lattice rule in \(t\) dimensions, whose points \(\mathbb{u}_i\) are defined by \( \mathbf{u}_i = (i \mathbf{a} \bmod m)/m \) where \(\mathbf{a} = (a_1,a_2,\dots,a_t) \in \mathbb{Z}_m^t\) is the generating vector; see Section 5.4 of the guide. More...
 
class  ReducerBB
 This class provides functions to find a shortest nonzero vector in the lattice using a BB algorithm as in [4], and to compute a Minkowski basis reduction as in [1]. More...
 
class  Weights
 Abstract class that defines an interface to specify Weights given to projections in figures of merit. More...
 
class  WeightsOrderDependent
 Defines order-dependent weights, for which the weight of a projection depends only on its order (cardinality). More...
 
class  WeightsPOD
 Defines product and order-dependent (POD) weights, for which the weight of a projection is the sum of a product weight and an order-dependent weight. More...
 
class  WeightsProduct
 Defines product weights, for which the weight of a projection is equal to the product of the weights of the individual coordinates. More...
 
class  WeightsProjectionDependent
 Defines projection-dependent weights, for which the weight for any given projection can be set individually by setWeight(). More...
 
class  WeightsUniform
 Specifies weights that are the same (usually 1) for all projections. More...
 

Typedefs

typedef double Weight
 A scalar weight type.
 

Enumerations

enum  NormType { SUPNORM , L1NORM , L2NORM , ZAREMBANORM }
 The available norm types to measure the length of vectors. More...
 
enum  OutputType { TERM , RES , TEX , GEN }
 Different choices of output formats. More...
 
enum  ProblemType {
  BASIS , DUAL , REDUCTION , SHORTEST ,
  MERIT
}
 Types of problems that LatticeTester can handle. More...
 
enum  PrecisionType { DOUBLE , QUADRUPLE , XDOUBLE , RR }
 This can be supersed by the Real type. More...
 
enum  PrimeType { PRIME , PROB_PRIME , COMPOSITE , UNKNOWN }
 Indicates whether an integer is prime, probably prime, composite or its status is unknown (or we do not care). More...
 
enum  CriterionType {
  LENGTH , SPECTRAL , BEYER , PALPHA ,
  BOUND_JS
}
 Merit criteria to measure the quality of generators or lattices. More...
 
enum  NormaType {
  BESTLAT , BESTUPBOUND , LAMINATED , ROGERS ,
  MINKL1 , MINKHLAW , NONE
}
 Different types of normalizations that can be used for shortest-vector lengths. More...
 
enum  CalcType { PAL , NORMPAL , BAL , SEEKPAL }
 Indicates which type of calculation is considered for the \(P_{\alpha}\) test. More...
 
enum  ReductionType {
  PAIR , LLL , BKZ , BB ,
  PAIRBB , LLLBB , BKZBB
}
 A list of all the possible lattice reductions implemented in LatticeTester. More...
 
enum  DecompTypeBB { CHOLESKY , TRIANGULAR }
 Two possible ways of obtaining a triangular matrix to compute the bounds in the BB algorithm. More...
 
enum  ProjConstructType { LLLPROJ , UPPERTRIPROJ }
 Two possible ways of computing the basis for a projection. More...
 
enum  MeritType { MERITM , MERITQ }
 Two different types of figures of merit. More...
 

Functions

template<typename Int , typename Real >
static long LLLConstruction0 (IntMat &gen, const double delta=0.9, long r=0, long c=0, RealVec *sqlen=0)
 This function takes a set of generating vectors of a lattice in matrix gen and finds a lattice basis by applying LLL reduction with the given value of delta, using the NTL implementation specified by prec.
 
template<typename Int , typename Real >
static void LLLBasisConstruction (IntMat &gen, const Int &m, const double delta=0.9, long r=0, long c=0, RealVec *sqlen=0)
 Similar to LLLConstruction0, except that this function adds implicitly the vectors \(m \mathbf{e}_i\) to the generating set, so it always returns a square matrix.
 
template<typename Int >
static void lowerTriangularBasis (IntMat &basis, IntMat &gen, const Int &m, long r=0, long c=0)
 Takes a set of generating vectors in the matrix gen and iteratively transforms it into a lower triangular lattice basis into the matrix basis.
 
template<typename Int >
static void upperTriangularBasis (IntMat &basis, IntMat &gen, const Int &m, long r=0, long c=0)
 Same as lowerTriangularBasis, except that the returned basis is upper triangular.
 
template<typename Int >
static void upperTriangularBasisOld96 (IntMat &basis, IntMat &gen, const Int &m, long r=0, long c=0)
 The old version from [3] and [6].
 
template<typename Int >
static void mDualLowerTriangular (IntMat &basisDual, const IntMat &basis, const Int &m, long dim=0)
 Takes a lower-triangular basis matrix basis and computes the m-dual basis basisDual.
 
template<typename Int >
static void mDualUpperTriangular (IntMat &basisDual, const IntMat &basis, const Int &m, long dim=0)
 Takes a upper triangular basis matrix basis and computes the m-dual basis basisDual.
 
template<typename Int >
static void mDualUpperTriangularOld96 (IntMat &basisDual, const IntMat &basis, const Int &m, long dim=0)
 This function does essentially the same thing as mDualUpperTriangular, but the algorithm is slightly different.
 
template<typename Int >
static void mDualBasis (IntMat &basisDual, const IntMat &basis, const Int &m)
 This function assumes that basis contains a basis of the primal lattice scaled by the factor m, not necessarily triangular, and it returns in basisDual the m-dual basis.
 
template<typename Int >
static void projectMatrix (IntMat &out, const IntMat &in, const Coordinates &proj, long r=0)
 This function overwrites the first r rows of matrix 'out' by a matrix formed by the first r rows of the c columns of matrix in that are specified by proj, where ‘c = size(proj)’ is the cardinality of the projection proj.
 
template<typename Int , typename Real >
static void projectionConstructionLLL (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, const double delta=0.9, long r=0, RealVec *sqlen=0)
 Constructs a basis for the projection proj of the lattice with basis inBasis, using LLLBasisConstruction, and returns it in projBasis.
 
template<typename Int >
static void projectionConstructionUpperTri (IntMat &projBasis, const IntMat &inBasis, IntMat &genTemp, const Coordinates &proj, const Int &m, long r=0)
 Same as projectionConstructionLLL, but the construction is made using upperTriangularBasis, so the returned basis is upper triangular.
 
template<typename Int >
static void projectionConstructionUpperTri (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, long r=0)
 
template<typename Int >
static void projectionConstruction (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, const ProjConstructType projType=LLLPROJ, const double delta=0.9)
 In this version, the construction method is passed as a parameter.
 
template<>
long LLLConstruction0 (NTL::Mat< long > &gen, const double delta, long r, long c, NTL::Vec< double > *sqlen)
 
template<>
long LLLConstruction0 (NTL::Mat< NTL::ZZ > &gen, const double delta, long r, long c, NTL::Vec< double > *sqlen)
 
template<>
long LLLConstruction0 (NTL::Mat< NTL::ZZ > &gen, const double delta, long r, long c, NTL::Vec< xdouble > *sqlen)
 
template<>
long LLLConstruction0 (NTL::Mat< NTL::ZZ > &gen, const double delta, long r, long c, NTL::Vec< quad_float > *sqlen)
 
template<>
long LLLConstruction0 (NTL::Mat< NTL::ZZ > &gen, const double delta, long r, long c, NTL::Vec< NTL::RR > *sqlen)
 
template<typename Int , typename Real >
void LLLBasisConstruction (IntMat &gen, const Int &m, double delta, long r, long c, RealVec *sqlen)
 Similar to LLLConstruction0, except that this function adds implicitly the vectors \(m \mathbf{e}_i\) to the generating set, so it always returns a square matrix.
 
template<typename Int >
void lowerTriangularBasis (IntMat &basis, IntMat &gen, const Int &m, long dim1, long dim2)
 Takes a set of generating vectors in the matrix gen and iteratively transforms it into a lower triangular lattice basis into the matrix basis.
 
template<typename Int >
void upperTriangularBasis (IntMat &basis, IntMat &gen, const Int &m, long dim1, long dim2)
 Same as lowerTriangularBasis, except that the returned basis is upper triangular.
 
template<typename Matr , typename Int >
void upperTriangularBasisOld96 (Matr &V, Matr &W, const Int &m, int64_t lin, int64_t col)
 
template<typename Int >
void mDualLowerTriangular (IntMat &B, const IntMat &A, const Int &m, long dim)
 Takes a lower-triangular basis matrix basis and computes the m-dual basis basisDual.
 
template<typename Int >
void mDualUpperTriangular (IntMat &B, const IntMat &A, const Int &m, long dim)
 Takes a upper triangular basis matrix basis and computes the m-dual basis basisDual.
 
template<typename Int >
void mDualUpperTriangularOld96 (IntMat &basisDual, const IntMat &basis, const Int &m, long dim)
 This function does essentially the same thing as mDualUpperTriangular, but the algorithm is slightly different.
 
template<>
void mDualUpperTriangularOld96 (NTL::Mat< NTL::ZZ > &basisDual, const NTL::Mat< NTL::ZZ > &basis, const NTL::ZZ &m, long dim)
 
template<typename Int >
void mDualBasis (NTL::Mat< Int > &basisDual, const NTL::Mat< Int > &basis, const Int &m)
 
template<>
void mDualBasis (NTL::Mat< NTL::ZZ > &basisDual, const NTL::Mat< NTL::ZZ > &basis, const NTL::ZZ &m)
 
template<typename Int >
void projectMatrix (IntMat &out, const IntMat &in, const Coordinates &proj, long r)
 This function overwrites the first r rows of matrix 'out' by a matrix formed by the first r rows of the c columns of matrix in that are specified by proj, where ‘c = size(proj)’ is the cardinality of the projection proj.
 
template<typename Int , typename Real >
void projectionConstructionLLL (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, const double delta, long r, RealVec *sqlen)
 Constructs a basis for the projection proj of the lattice with basis inBasis, using LLLBasisConstruction, and returns it in projBasis.
 
template<typename Int >
void projectionConstructionUpperTri (IntMat &projBasis, const IntMat &inBasis, IntMat &genTemp, const Coordinates &proj, const Int &m, long r)
 Same as projectionConstructionLLL, but the construction is made using upperTriangularBasis, so the returned basis is upper triangular.
 
template<typename Int >
void projectionConstructionUpperTri (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, long r)
 
template<typename Int >
void projectionConstruction (IntMat &projBasis, const IntMat &inBasis, const Coordinates &proj, const Int &m, const ProjConstructType projType, const double delta)
 In this version, the construction method is passed as a parameter.
 
static std::string toStringNorm (NormType norm)
 The following are functions for printing the enum constants in this module.
 
static std::string toStringOutput (OutputType out)
 
static std::string toStringProblem (ProblemType prob)
 
static std::string toStringPrecision (PrecisionType precision)
 
static std::string toStringPrime (PrimeType prim)
 
static std::string toStringCriterion (CriterionType criter)
 
static std::string toStringNorma (NormaType norma)
 
static std::string toStringCalc (CalcType calc)
 
static std::string toStringReduction (ReductionType reduct)
 
static std::string toStringDecomp (DecompTypeBB decomp)
 
static std::string toStringProjConstruct (ProjConstructType proj)
 
static std::string toStringMeritType (MeritType merit)
 
template<typename Int , typename Real >
static void redLLL (IntMat &basis, double delta=0.99999, long dim=0, RealVec *sqlen=0)
 This function uses the NTL implementation of the LLL reduction algorithm with factor delta, presented in [9] (see also [7]).
 
template<typename Int >
static void redLLLExact (IntMat &basis, double delta=0.99999)
 This static function implements an exact algorithm from NTL to perform the original LLL reduction.
 
template<typename Int , typename Real >
static void redBKZ (IntMat &basis, double delta=0.99999, int64_t blocksize=10, long prune=0, long dim=0, RealVec *sqlen=0)
 This calls the NTL implementation of the floating point version of the BKZ reduction algorithm presented in [9], with reduction factor delta, block size blocksize, pruning parameter ‘prune’; see [7].
 
template<typename Int , typename Real >
void redLLL (IntMat &basis, double delta, long dim, RealVec *sqlen)
 This function uses the NTL implementation of the LLL reduction algorithm with factor delta, presented in [9] (see also [7]).
 
template<>
void redLLL (NTL::Mat< int64_t > &basis, double delta, long dim, NTL::Vec< double > *sqlen)
 
template<>
void redLLL (NTL::Mat< NTL::ZZ > &basis, double delta, long dim, NTL::Vec< double > *sqlen)
 
template<>
void redLLL (NTL::Mat< NTL::ZZ > &basis, double delta, long dim, NTL::Vec< xdouble > *sqlen)
 
template<>
void redLLL (NTL::Mat< NTL::ZZ > &basis, double delta, long dim, NTL::Vec< quad_float > *sqlen)
 
template<>
void redLLL (NTL::Mat< NTL::ZZ > &basis, double delta, long dim, NTL::Vec< NTL::RR > *sqlen)
 
template<typename Int >
void redLLLExact (IntMat &basis, double delta)
 This static function implements an exact algorithm from NTL to perform the original LLL reduction.
 
template<>
void redLLLExact (NTL::Mat< NTL::ZZ > &basis, double delta)
 
template<typename Int , typename Real >
void redBKZ (IntMat &basis, double delta, long blocksize, long prune, long dim, RealVec *sqlen)
 
template<>
void redBKZ (NTL::Mat< int64_t > &basis, double delta, long blocksize, long prune, long dim, NTL::Vec< double > *sqlen)
 
template<>
void redBKZ (NTL::Mat< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::Vec< double > *sqlen)
 
template<>
void redBKZ (NTL::Mat< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::Vec< xdouble > *sqlen)
 
template<>
void redBKZ (NTL::Mat< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::Vec< quad_float > *sqlen)
 
template<>
void redBKZ (NTL::Mat< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::Vec< NTL::RR > *sqlen)
 
template<typename T >
void swap9 (T &x, T &y)
 Takes references to two variables of a generic type and swaps their content.
 
template<typename Matr1 , typename Matr2 >
void copyMatrixToMat (Matr1 &A, Matr2 &B)
 
template<typename Int >
void printBase (IntMat bas_mat)
 
template<typename Int >
void printBase2 (IntMat bas_mat)
 
template<typename Int >
void copy (IntMat &b1, IntMat &b2)
 
template<typename Int >
void copy (IntMat &b1, IntMat &b2, int64_t r, int64_t c)
 
int64_t getWidth (clock_t time[], int64_t dim, std::string message, clock_t totals[], int64_t ind)
 
static std::ostream & operator<< (std::ostream &os, const Coordinates &coords)
 Formats the coordinate set coords and outputs it to os.
 
static std::istream & operator>> (std::istream &is, Coordinates &coords)
 Reads a formatted coordinate set from is.
 
std::ostream & operator<< (std::ostream &os, const Weights &o)
 Identifies the type of weights, formats them and outputs them on os.
 
std::istream & operator>> (std::istream &is, WeightsProjectionDependent &weights)
 Reads formatted projection-dependent weights into the object weights.
 
Random numbers

All the functions of this module use LFSR258 as an underlying source for pseudo-random numbers.

A free (as in freedom) implementation of this generator can be found at http://simul.iro.umontreal.ca/ as well as the article presenting it. All the functions generating some sort of random number will advance an integer version of LFSR258 by one state and output a transformation of the state to give a double, an int64_t or bits.

double RandU01 ()
 Returns a random number in \([0, 1)\).
 
int64_t RandInt (int64_t i, int64_t j)
 Return a uniform pseudo-random integer in \([i, j]\).
 
std::uint64_t RandBits (int64_t s)
 Returns the first s pseudo-random bits of the underlying RNG in the form of a s-bit integer.
 
void SetSeed (std::uint64_t seed)
 Sets the seed of the generator.
 
Mathematical functions

These are complete reimplementation of certain mathematical functions, or wrappers for standard C/C++ functions.

double mysqrt (double x)
 Returns \(\sqrt{x}\) for \(x\ge0\), and \(-1\) for \(x < 0\).
 
template<typename T >
double Lg (const T &x)
 Returns the logarithm of \(x\) in base 2.
 
double Lg (std::int64_t x)
 Returns the logarithm of \(x\) in base 2.
 
template<typename Scal >
Scal abs (Scal x)
 Returns the absolute value of \(x\).
 
template<typename T >
std::int64_t sign (const T &x)
 Returns the sign of x.
 
template<typename Real >
Real Round (Real x)
 Returns the value of x rounded to the NEAREST integer value.
 
std::int64_t Factorial (int64_t t)
 Calculates \(t!\), the factorial of \(t\) and returns it as an std::int64_t.
 
Division and modular arithmetic
Remarks
Richard: Pour certaines fonctions, les résultats sont mis dans les premiers arguments de la fonction pour être compatible avec NTL; pour d’autres, ils sont mis dans les derniers arguments pour être compatible avec notre ancienne version de LatMRG en Modula-2. Plutôt détestable. Je crois qu’il faudra un jour réarranger les arguments des fonctions pour qu’elles suivent toutes la même convention que NTL.

This module offers function to perform division and find remainders in a standard way. These functions are usefull in the case where one wants to do divisions or find remainders of operations with negative operands. The reason is that NTL and primitive types do not use the same logic when doing calculations on negative numbers.

Basically, NTL will always floor a division and C++ will always truncate a division (which effectively means the floor function is replaced by a roof function if the answer is a negative number). When calculating the remainder of x/y, both apply the same logic but get a different result because they do not do the same division. In both representations, we have that

\[ y\cdot(x/y) + x%y = x. \]

It turns out that, with negative values, NTL will return an integer with the same sign as y where C++ will return an integer of opposite sign (but both will return the same number modulo y).

template<typename Int >
void Quotient (const Int &a, const Int &b, Int &q)
 Computes a/b, truncates the fractional part and puts the result in q.
 
template<typename Real >
void Modulo (const Real &a, const Real &b, Real &r)
 Computes the remainder of a/b and stores its positive equivalent mod b in r.
 
template<typename Int >
void ModuloTowardZero (const Int &a, const Int &b, Int &r)
 Computes the element r = a mod b that is closest to 0.
 
void ModuloTowardZero (const std::int64_t &a, const std::int64_t &b, std::int64_t &r)
 
template<typename Real >
void Divide (Real &q, Real &r, const Real &a, const Real &b)
 Computes the quotient \(q = a/b\) and remainder \(r = a \bmod b\).
 
template<typename Real >
void DivideRound (const Real &a, const Real &b, Real &q)
 Computes \(a/b\), rounds the result to the nearest integer and returns the result in \(q\).
 
std::int64_t gcd (std::int64_t a, std::int64_t b)
 Returns the value of the greatest common divisor of \(a\) and \(b\) by using Stein's binary GCD algorithm.
 
template<typename Int >
void Euclide (const Int &A, const Int &B, Int &C, Int &D, Int &E, Int &F, Int &G)
 This method computes the greater common divisor of A and B with Euclid's algorithm.
 
template<typename Int >
void Euclide (const Int &A, const Int &B, Int &C, Int &D, Int &G)
 
template<typename Int >
void TransposeMatrix (IntMat &mat, IntMat &mat2)
 Takes a set of generating vectors in the matrix mat and iteratively transforms it into an upper triangular lattice basis into the matrix mat2.
 
Vectors

These are utilities to manipulate vectors ranging from instantiation to scalar product.

******** MOST OF THIS CAN BE REMOVED. WE CAN JUST USE NTL DIRECTLY. *********

template<typename Real >
void CreateVect (Real *&A, int64_t d)
 Allocates memory to A as an array of Real of dimension d and initializes its elements to 0.
 
template<typename Vect >
void CreateVect (Vect &A, int64_t d)
 Creates the vector A of dimensions d+1 and initializes its elements to 0.
 
template<typename Real >
void DeleteVect (Real *&A)
 Frees the memory used by the vector A.
 
template<typename Vect >
void DeleteVect (Vect &A)
 Frees the memory used by the vector A, destroying all the elements it contains.
 
template<typename Real >
void SetZero (Real *A, int64_t d)
 Sets the first d of A to 0.
 
template<typename Vect >
void SetZero (Vect &A, int64_t d)
 Sets the first d components of A to 0.
 
template<typename Real >
void SetValue (Real *A, int64_t d, const Real &x)
 Sets the first d components of A to 0.
 
template<typename Vect >
std::string toString (const Vect &A, int64_t c, int64_t d, const char *sep=" ")
 Returns a string containing A[c] to A[d-1] formated as [A[c]sep...sepA[d-1]].
 
template<typename Vect >
std::string toString (const Vect &A, int64_t d)
 Returns a string containing the first d components of the vector A as a string.
 
template<typename Int , typename Real >
void ProdScal (const IntVec &A, const IntVec &B, int64_t n, Real &D)
 Computes the scalar product of vectors A and B truncated to their n first components, then puts the result in D.
 
template<typename Int >
void Invert (const IntVec &A, IntVec &B, int64_t n)
 Takes an input vector A of dimension n+1 and fill the vector B with the values [-A[n] -A[n-1] ... -A[1][1].
 
template<typename Int , typename Real >
void CalcNorm (const IntVec &V, int64_t n, Real &S, NormType norm)
 Computes the norm norm of vector V trunctated to its n first components, and puts the result in S.
 
template<typename Vect >
void CopyPartVec (Vect &toVec, const Vect &fromVec, int64_t c)
 Copies the first c components of vector fromVec into vector toVec.
 
template<typename Matr >
void CopyPartMat (Matr &toMat, const Matr &fromMat, int64_t r, int64_t c)
 Copies the first r rows and c columns of matrix fromMat into matrix toMat.
 
template<typename Vect >
void CopyVect (Vect &A, const Vect &B, int64_t n)
 Copies the first n components of vector B into vector A.
 
template<typename Vect , typename Scal >
void ModifVect (Vect &A, const Vect &B, Scal x, int64_t n)
 Adds the first n components of vector B multiplied by x to first n components of vector A.
 
template<typename Vect , typename Scal , typename Int >
void ModifVectModulo (Vect &A, const Vect &B, Scal x, Int m, int64_t n)
 Adds the first n components of vector B multiplied by x to first n components of vector A, modulo m.
 
template<typename Vect >
void ChangeSign (Vect &A, int64_t n)
 Changes the sign (multiplies by -1) the first n components of vector A.
 
std::int64_t GCD2vect (std::vector< std::int64_t > V, int64_t k, int64_t n)
 Computes the greatest common divisor of V[k],...,V[n-1].
 
Matrices
template<typename Real >
void CreateMatr (Real **&A, int64_t d)
 Allocates memory to a square matrix A of dimensions \(d \times d\) and initializes its elements to 0.
 
template<typename Real >
void CreateMatr (Real **&A, int64_t line, int64_t col)
 Allocates memory for the matrix A of dimensions \(\text{line} \times \text{col}\) and initializes its elements to 0.
 
template<typename Int >
void CreateMatr (IntMat &A, int64_t d)
 Resizes the matrix A to a square matrix of dimensions d*d and re-initializes its elements to 0.
 
template<typename Int >
void CreateMatr (IntMat &A, int64_t line, int64_t col)
 Resizes the matrix A to a matrix of dimensions \(line \times col\) and re-initializes its elements to 0.
 
template<typename Real >
void DeleteMatr (Real **&A, int64_t d)
 Frees the memory used by the \(d \times d\) matrix A.
 
template<typename Real >
void DeleteMatr (Real **&A, int64_t line, int64_t col)
 Frees the memory used by the matrix A of dimension \(\text{line} \times \text{col}\).
 
template<typename Int >
void DeleteMatr (IntMat &A)
 Calls the clear() method on A.
 
template<typename Matr >
void CopyMatr (Matr &A, const Matr &B, int64_t n)
 Copies the \(n \times n\) submatrix of the first lines and columns of B into matrix A.
 
template<typename Matr >
void CopyMatr (Matr &A, const Matr &B, int64_t line, int64_t col)
 Copies the \(\text{line} \times col\) submatrix of the first lines and columns of B into matrix A.
 
template<typename MatT >
std::string toStr (const MatT &mat, int64_t d1, int64_t d2, int64_t prec=2)
 Returns a string that is a representation of mat.
 
template<typename Int >
void ProductDiagonal (const NTL::Mat< Int > &A, long dim, Int &prod)
 Returns the product of the diagonal elements of the matrix A, which is assumed to be square dim x dim.
 
template<typename Int >
bool CheckTriangular (const NTL::Mat< Int > &A, long dim, const Int m)
 Checks that the upper \(\text{dim} \times \text{dim}\) submatrix of A is triangular modulo m.
 
template<typename Int >
bool checkInverseModm (const NTL::Mat< Int > &A, const NTL::Mat< Int > &B, const Int m)
 Checks if A and B are m-dual to each other.
 
template<typename Matr , typename Int >
void calcDual (const Matr &A, Matr &B, int64_t d, const Int &m)
 Takes an upper triangular basis A and computes an m-dual lattice basis to this matrix.
 
Debugging functions
void myExit (std::string msg)
 Simple error exit function, prints msg on exit.
 
Printing functions and operators
template<class K , class T , class C , class A >
std::ostream & operator<< (std::ostream &out, const std::map< K, T, C, A > &x)
 Streaming operator for maps.
 
template<class T1 , class T2 >
std::ostream & operator<< (std::ostream &out, const std::pair< T1, T2 > &x)
 Streaming operator for vectors.
 
template<class T , class A >
std::ostream & operator<< (std::ostream &out, const std::vector< T, A > &x)
 Streaming operator for vectors.
 
template<class K , class C , class A >
std::ostream & operator<< (std::ostream &out, const std::set< K, C, A > &x)
 Streaming operator for sets.
 

Variables

const double MAX_LONG_DOUBLE = 9007199254740992.0
 Maximum integer that can be represented exactly as a double: \(2^{53}\).
 
const std::int64_t TWO_EXP []
 Table of powers of 2: TWO_EXP[ \(i\)] \(= 2^i\), \(i=0, 1, …, 63\).
 

Detailed Description

LatticeTester namespace.

Typedef Documentation

◆ Weight

typedef double LatticeTester::Weight

A scalar weight type.

Enumeration Type Documentation

◆ NormType

The available norm types to measure the length of vectors.

For \(X = (x_1,…,x_t)\):
SUPNORM corresponds to \(\Vert X\Vert= \max(|x_1|,…,|x_t|)\).
L1NORM corresponds to \(\Vert X\Vert= |x_1|+\cdots+|x_t|\).
L2NORM corresponds to \(\Vert X\Vert= (x_1^2+\cdots+x_t^2)^{1/2}\).
ZAREMBANORM corresponds to \(\Vert X\Vert= \max(1, |x_1|)\cdots\max(1, |x_t|)\).

Enumerator
SUPNORM 
L1NORM 
L2NORM 
ZAREMBANORM 

◆ OutputType

Different choices of output formats.

TERM: the results will appear only on the terminal screen.
RES: the results will be in plain text and sent to a .res file.
TEX: the results will be in a LaTeX file with extension .tex.
GEN: a list of retained generators will be sent to a file with extension .gen, in a specific format so this list can be read again for further analysis.

Enumerator
TERM 
RES 
TEX 
GEN 

◆ ProblemType

Types of problems that LatticeTester can handle.

Not sure if we still need this. *******

Enumerator
BASIS 
DUAL 
REDUCTION 
SHORTEST 
MERIT 

◆ PrecisionType

This can be supersed by the Real type.

Types of precision that the NTL can use for real numbers: DOUBLE – double QUADRUPLE – quad_float (quasi quadruple precision) this is useful when roundoff errors can cause problems XDOUBLE – xdouble (extended exponent doubles) this is useful when numbers get too big RR – RR (arbitrary precision floating point). The choice DOUBLE is usually the fastest, but may be prone to roundoff errors and/or overflow. See https://github.com/u-u-h/NTL/blob/master/doc/LLL.txt.

Enumerator
DOUBLE 
QUADRUPLE 
XDOUBLE 
RR 

◆ PrimeType

Indicates whether an integer is prime, probably prime, composite or its status is unknown (or we do not care).

Enumerator
PRIME 
PROB_PRIME 
COMPOSITE 
UNKNOWN 

◆ CriterionType

Merit criteria to measure the quality of generators or lattices.

TO DO: this list is not very clear. Maybe outdated. ****************

LENGTH: Only using the length of the shortest vector as a criterion. SPECTRAL: figure of merit \(S_T\) based on the spectral test.
BEYER: figure of merit is the Beyer quotient \(Q_T\).
PALPHA: figure of merit based on \(P_{\alpha}\).
BOUND_JS: figure of merit based on the Joe-Sinescu bound [10].
???

Enumerator
LENGTH 
SPECTRAL 
BEYER 
PALPHA 
BOUND_JS 

◆ NormaType

Different types of normalizations that can be used for shortest-vector lengths.

Corresponds to different ways of approximating the Hermite constants gamma_t.

BESTLAT: the value used for \(d_t^*\) corresponds to the best lattice.
BESTUPBOUND: the value used for \(d_t^*\) corresponds to the best bound known to us.
LAMINATED: the value used for \(d_t^*\) corresponds to the best laminated lattice.
ROGERS: the value for \(d_t^*\) is obtained from Rogers’ bound on the density of sphere packing.
MINKL1: the value for \(d_t^*\) is obtained from the theoretical bounds on the length of the shortest nonzero vector in the lattice using the \({\mathcal{L}}_1\) norm.
MINKHLAW: the value for \(d_t^*\) is obtained from Minkowski’ theoretical bounds on the length of the shortest nonzero vector in the lattice using the \({\mathcal{L}}_2\) norm.
NONE: no normalization will be used.

Enumerator
BESTLAT 
BESTUPBOUND 
LAMINATED 
ROGERS 
MINKL1 
MINKHLAW 
NONE 

◆ CalcType

Indicates which type of calculation is considered for the \(P_{\alpha}\) test.

Is this used anywhere? ************

PAL is for the \(P_{\alpha}\) test.
BAL is for the bound on the \(P_{\alpha}\) test.
NORMPAL is for the \(P_{\alpha}\) test PAL, with the result normalized over the BAL bound.
SEEKPAL is for the \(P_{\alpha}\) seek, which searches for good values of the multiplier.

Enumerator
PAL 
NORMPAL 
BAL 
SEEKPAL 

◆ ReductionType

A list of all the possible lattice reductions implemented in LatticeTester.

PAIR: Pairwise reductions only. LLL: LLL reduction only. BKZ: block Korkine-Zolotarev reduction only. BB: direct shortest vector search with BB (no pre-red.). PAIRBB: Pairwise reduction followed by BB. LLLBB: LLL followed by BB. BKZBB: BKZ followed by BB.

Enumerator
PAIR 
LLL 
BKZ 
BB 
PAIRBB 
LLLBB 
BKZBB 

◆ DecompTypeBB

Two possible ways of obtaining a triangular matrix to compute the bounds in the BB algorithm.

CHOLESKY: use a lower-triangular matrix obtained as the Cholesky decomposition of the matrix of scalar products. TRIANGULAR: use a lower-triangular basis

Enumerator
CHOLESKY 
TRIANGULAR 

◆ ProjConstructType

Two possible ways of computing the basis for a projection.

LLLPROJ: uses LLL reduction. UPPERTRIPROJ: use an upper-triangular basis construction.

Enumerator
LLLPROJ 
UPPERTRIPROJ 

◆ MeritType

Two different types of figures of merit.

MERITM: based on shortest vector. MERITQ: based on Beyer quotient.

Enumerator
MERITM 
MERITQ 

Function Documentation

◆ LLLConstruction0() [1/6]

template<typename Int , typename Real >
static long LatticeTester::LLLConstruction0 ( IntMat & gen,
const double delta = 0.9,
long r = 0,
long c = 0,
RealVec * sqlen = 0 )
static

This function takes a set of generating vectors of a lattice in matrix gen and finds a lattice basis by applying LLL reduction with the given value of delta, using the NTL implementation specified by prec.

The basis is returned in the first rows of gen, in a number of rows equal to its rank. See the file EnumTypes and the documentation of LLL in NTL for the meaning and choices for prec. The default value of delta is not very close to 1 because the goal here is just to compute a basis. It can be taken much closer to 1 if we really prefer a highly-reduced basis. The optional parameters r and c indicate the numbers of rows and columns of the matrix gen that are actually used. When they are 0 (the default values), then these numbers are taken to be the dimensions of the matrix gen. When sqlen is not 0, the square lengths of the basis vectors are returned in this array, exactly as in the LLL_lt.h module. The function returns the dimension (number of rows) of the newly computed basis, which may differ from the number of rows of the gen object. The latter is never resized.

This function does not assume that all vectors \(m e_i\) belong to the lattice, so it may return a basis matrix that has fewer rows than columns! To make sure that these vectors belong to the lattice, we can add them explicitly beforehand to the set of generating vectors, or call the next function.

◆ LLLBasisConstruction() [1/2]

template<typename Int , typename Real >
static void LatticeTester::LLLBasisConstruction ( IntMat & gen,
const Int & m,
const double delta = 0.9,
long r = 0,
long c = 0,
RealVec * sqlen = 0 )
static

Similar to LLLConstruction0, except that this function adds implicitly the vectors \(m \mathbf{e}_i\) to the generating set, so it always returns a square matrix.

The matrix gen is not resized by this function, so it can remain larger than the lattice dimension.

◆ lowerTriangularBasis() [1/2]

template<typename Int >
static void LatticeTester::lowerTriangularBasis ( IntMat & basis,
IntMat & gen,
const Int & m,
long r = 0,
long c = 0 )
static

Takes a set of generating vectors in the matrix gen and iteratively transforms it into a lower triangular lattice basis into the matrix basis.

This lattice is assumed to contain all the vectors of the form \(m \mathbf{e}_j\), so these vectors are added implicitly to the generating set. Apart from that, all the entries of gen given as input are assumed to be reduced modulo the scaling factor m and all the computations are done modulo m. The matrix basis is assumed to be large enough to contain the new basis, whose dimension should be the number of columns taken from gen. After the execution, gen will contain irrelevant information (garbage) and basis will contain an upper triangular basis. The algorithm is explained in the lattice tester guide. Important: gen and basis must be different IntMat objects.

When r and/or c are strictly positive, they specify the numbers of rows and columns of gen that are actually used, as in LLLConstruction0. The matrix gen is never resized. The new basis should be c x c.

◆ upperTriangularBasis() [1/2]

template<typename Int >
static void LatticeTester::upperTriangularBasis ( IntMat & basis,
IntMat & gen,
const Int & m,
long r = 0,
long c = 0 )
static

Same as lowerTriangularBasis, except that the returned basis is upper triangular.

◆ upperTriangularBasisOld96() [1/2]

template<typename Int >
static void LatticeTester::upperTriangularBasisOld96 ( IntMat & basis,
IntMat & gen,
const Int & m,
long r = 0,
long c = 0 )
static

The old version from [3] and [6].

◆ mDualLowerTriangular() [1/2]

template<typename Int >
static void LatticeTester::mDualLowerTriangular ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )
static

Takes a lower-triangular basis matrix basis and computes the m-dual basis basisDual.

The function assumes that each coefficient on the diagonal of basis is nonzero and divides m. That is, the basis matrix must be square and m-invertible. Since the basis is lower triangular, its m-dual will be upper triangular. When dim > 0, it must give the number of rows and columns of the matrix basis that is actually used. Otherwise (by default) the function uses basis.numCols().

◆ mDualUpperTriangular() [1/2]

template<typename Int >
static void LatticeTester::mDualUpperTriangular ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )
static

Takes a upper triangular basis matrix basis and computes the m-dual basis basisDual.

This function is the equivalent of mDualLowerTriangular for upper-triangular matrices.

◆ mDualUpperTriangularOld96() [1/3]

template<typename Int >
static void LatticeTester::mDualUpperTriangularOld96 ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )
static

This function does essentially the same thing as mDualUpperTriangular, but the algorithm is slightly different.

It uses the method described in [3].

◆ mDualBasis() [1/3]

template<typename Int >
static void LatticeTester::mDualBasis ( IntMat & basisDual,
const IntMat & basis,
const Int & m )
static

This function assumes that basis contains a basis of the primal lattice scaled by the factor m, not necessarily triangular, and it returns in basisDual the m-dual basis.

It uses matrix inversion and is rather slow. It is currently implemented only for Int = ZZ and it also assumes that the dimensions of the two IntMat objects is exactly the same as the dimensions of the lattices. The reason for this is that we use an NTL function that works only under these conditions.

◆ projectMatrix() [1/2]

template<typename Int >
static void LatticeTester::projectMatrix ( IntMat & out,
const IntMat & in,
const Coordinates & proj,
long r = 0 )
static

This function overwrites the first r rows of matrix 'out' by a matrix formed by the first r rows of the c columns of matrix in that are specified by proj, where ‘c = size(proj)’ is the cardinality of the projection proj.

Only the first c columns of the first r rows are overwritten; the other entries are left unchanged. If r = 0, then all the rows of matrix in are taken. When in is a basis, r will usually be the dimension of that basis and c will be smaller. The dimensions of the matrices in and out are always left unchanged. The dimensions of out are assumed to be large enough. After the call, the matrix out will then contain a set of generating vectors for the projection proj. The matrices in and out must be different IntMat objects, otherwise the program halts with an error message.

◆ projectionConstructionLLL() [1/2]

template<typename Int , typename Real >
static void LatticeTester::projectionConstructionLLL ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
const double delta = 0.9,
long r = 0,
RealVec * sqlen = 0 )
static

Constructs a basis for the projection proj of the lattice with basis inBasis, using LLLBasisConstruction, and returns it in projBasis.

This returned basis is not triangular in general. Its dimension will be the number of coordinates in proj. The matrix projBasis must have enough columns to hold it and at least as many rows as the number of rows that we use from inBasis. When r > 0, only the first r rows of the matrix inBasis are used, otherwise we use all rows of that matrix. This r should normally be the true dimension dim of that basis, which is often smaller than the size maxDim of the IntMat object that contains the basis. The square Euclidean lengths of the basis vectors are returned in the array sqlen when the latter is given.

◆ projectionConstructionUpperTri() [1/4]

template<typename Int >
static void LatticeTester::projectionConstructionUpperTri ( IntMat & projBasis,
const IntMat & inBasis,
IntMat & genTemp,
const Coordinates & proj,
const Int & m,
long r = 0 )
static

Same as projectionConstructionLLL, but the construction is made using upperTriangularBasis, so the returned basis is upper triangular.

When r > 0, only the first r rows of the matrix inBasis are actually used, otherwise we use all rows of that matrix. In the first version, we pass a matrix genTemp that will be used to store the generating vectors of the projection before making the triangularization. The two matrices projBasis and genTemp must have enough columns to hold the projection and at least as many rows as the number of rows that we use from inBasis. We pass genTemp as a parameter to avoid the internal creation of a new matrix each time, in case we call this function several times. Its contents will be modified. In the second version, this matrix is not passed and a temporary one is created internally, which may add a bit of overhead.

◆ projectionConstructionUpperTri() [2/4]

template<typename Int >
static void LatticeTester::projectionConstructionUpperTri ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
long r = 0 )
static

◆ projectionConstruction() [1/2]

template<typename Int >
static void LatticeTester::projectionConstruction ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
const ProjConstructType projType = LLLPROJ,
const double delta = 0.9 )
static

In this version, the construction method is passed as a parameter.

The default is LLL. In the triangular case, a temporary matrix is created internally.

◆ LLLConstruction0() [2/6]

template<>
long LatticeTester::LLLConstruction0 ( NTL::Mat< long > & gen,
const double delta,
long r,
long c,
NTL::Vec< double > * sqlen )

◆ LLLConstruction0() [3/6]

template<>
long LatticeTester::LLLConstruction0 ( NTL::Mat< NTL::ZZ > & gen,
const double delta,
long r,
long c,
NTL::Vec< double > * sqlen )

◆ LLLConstruction0() [4/6]

template<>
long LatticeTester::LLLConstruction0 ( NTL::Mat< NTL::ZZ > & gen,
const double delta,
long r,
long c,
NTL::Vec< xdouble > * sqlen )

◆ LLLConstruction0() [5/6]

template<>
long LatticeTester::LLLConstruction0 ( NTL::Mat< NTL::ZZ > & gen,
const double delta,
long r,
long c,
NTL::Vec< quad_float > * sqlen )

◆ LLLConstruction0() [6/6]

template<>
long LatticeTester::LLLConstruction0 ( NTL::Mat< NTL::ZZ > & gen,
const double delta,
long r,
long c,
NTL::Vec< NTL::RR > * sqlen )

◆ LLLBasisConstruction() [2/2]

template<typename Int , typename Real >
void LatticeTester::LLLBasisConstruction ( IntMat & gen,
const Int & m,
const double delta = 0.9,
long r = 0,
long c = 0,
RealVec * sqlen = 0 )

Similar to LLLConstruction0, except that this function adds implicitly the vectors \(m \mathbf{e}_i\) to the generating set, so it always returns a square matrix.

The matrix gen is not resized by this function, so it can remain larger than the lattice dimension.

◆ lowerTriangularBasis() [2/2]

template<typename Int >
void LatticeTester::lowerTriangularBasis ( IntMat & basis,
IntMat & gen,
const Int & m,
long r = 0,
long c = 0 )

Takes a set of generating vectors in the matrix gen and iteratively transforms it into a lower triangular lattice basis into the matrix basis.

This lattice is assumed to contain all the vectors of the form \(m \mathbf{e}_j\), so these vectors are added implicitly to the generating set. Apart from that, all the entries of gen given as input are assumed to be reduced modulo the scaling factor m and all the computations are done modulo m. The matrix basis is assumed to be large enough to contain the new basis, whose dimension should be the number of columns taken from gen. After the execution, gen will contain irrelevant information (garbage) and basis will contain an upper triangular basis. The algorithm is explained in the lattice tester guide. Important: gen and basis must be different IntMat objects.

When r and/or c are strictly positive, they specify the numbers of rows and columns of gen that are actually used, as in LLLConstruction0. The matrix gen is never resized. The new basis should be c x c.

◆ upperTriangularBasis() [2/2]

template<typename Int >
void LatticeTester::upperTriangularBasis ( IntMat & basis,
IntMat & gen,
const Int & m,
long dim1,
long dim2 )

Same as lowerTriangularBasis, except that the returned basis is upper triangular.

◆ upperTriangularBasisOld96() [2/2]

template<typename Matr , typename Int >
void LatticeTester::upperTriangularBasisOld96 ( Matr & V,
Matr & W,
const Int & m,
int64_t lin,
int64_t col )

◆ mDualLowerTriangular() [2/2]

template<typename Int >
void LatticeTester::mDualLowerTriangular ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )

Takes a lower-triangular basis matrix basis and computes the m-dual basis basisDual.

The function assumes that each coefficient on the diagonal of basis is nonzero and divides m. That is, the basis matrix must be square and m-invertible. Since the basis is lower triangular, its m-dual will be upper triangular. When dim > 0, it must give the number of rows and columns of the matrix basis that is actually used. Otherwise (by default) the function uses basis.numCols().

◆ mDualUpperTriangular() [2/2]

template<typename Int >
void LatticeTester::mDualUpperTriangular ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )

Takes a upper triangular basis matrix basis and computes the m-dual basis basisDual.

This function is the equivalent of mDualLowerTriangular for upper-triangular matrices.

◆ mDualUpperTriangularOld96() [2/3]

template<typename Int >
void LatticeTester::mDualUpperTriangularOld96 ( IntMat & basisDual,
const IntMat & basis,
const Int & m,
long dim = 0 )

This function does essentially the same thing as mDualUpperTriangular, but the algorithm is slightly different.

It uses the method described in [3].

◆ mDualUpperTriangularOld96() [3/3]

template<>
void LatticeTester::mDualUpperTriangularOld96 ( NTL::Mat< NTL::ZZ > & basisDual,
const NTL::Mat< NTL::ZZ > & basis,
const NTL::ZZ & m,
long dim )

◆ mDualBasis() [2/3]

template<typename Int >
void LatticeTester::mDualBasis ( NTL::Mat< Int > & basisDual,
const NTL::Mat< Int > & basis,
const Int & m )

◆ mDualBasis() [3/3]

template<>
void LatticeTester::mDualBasis ( NTL::Mat< NTL::ZZ > & basisDual,
const NTL::Mat< NTL::ZZ > & basis,
const NTL::ZZ & m )

◆ projectMatrix() [2/2]

template<typename Int >
void LatticeTester::projectMatrix ( IntMat & out,
const IntMat & in,
const Coordinates & proj,
long r = 0 )

This function overwrites the first r rows of matrix 'out' by a matrix formed by the first r rows of the c columns of matrix in that are specified by proj, where ‘c = size(proj)’ is the cardinality of the projection proj.

Only the first c columns of the first r rows are overwritten; the other entries are left unchanged. If r = 0, then all the rows of matrix in are taken. When in is a basis, r will usually be the dimension of that basis and c will be smaller. The dimensions of the matrices in and out are always left unchanged. The dimensions of out are assumed to be large enough. After the call, the matrix out will then contain a set of generating vectors for the projection proj. The matrices in and out must be different IntMat objects, otherwise the program halts with an error message.

◆ projectionConstructionLLL() [2/2]

template<typename Int , typename Real >
void LatticeTester::projectionConstructionLLL ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
const double delta = 0.9,
long r = 0,
RealVec * sqlen = 0 )

Constructs a basis for the projection proj of the lattice with basis inBasis, using LLLBasisConstruction, and returns it in projBasis.

This returned basis is not triangular in general. Its dimension will be the number of coordinates in proj. The matrix projBasis must have enough columns to hold it and at least as many rows as the number of rows that we use from inBasis. When r > 0, only the first r rows of the matrix inBasis are used, otherwise we use all rows of that matrix. This r should normally be the true dimension dim of that basis, which is often smaller than the size maxDim of the IntMat object that contains the basis. The square Euclidean lengths of the basis vectors are returned in the array sqlen when the latter is given.

◆ projectionConstructionUpperTri() [3/4]

template<typename Int >
void LatticeTester::projectionConstructionUpperTri ( IntMat & projBasis,
const IntMat & inBasis,
IntMat & genTemp,
const Coordinates & proj,
const Int & m,
long r = 0 )

Same as projectionConstructionLLL, but the construction is made using upperTriangularBasis, so the returned basis is upper triangular.

When r > 0, only the first r rows of the matrix inBasis are actually used, otherwise we use all rows of that matrix. In the first version, we pass a matrix genTemp that will be used to store the generating vectors of the projection before making the triangularization. The two matrices projBasis and genTemp must have enough columns to hold the projection and at least as many rows as the number of rows that we use from inBasis. We pass genTemp as a parameter to avoid the internal creation of a new matrix each time, in case we call this function several times. Its contents will be modified. In the second version, this matrix is not passed and a temporary one is created internally, which may add a bit of overhead.

◆ projectionConstructionUpperTri() [4/4]

template<typename Int >
void LatticeTester::projectionConstructionUpperTri ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
long r )

◆ projectionConstruction() [2/2]

template<typename Int >
void LatticeTester::projectionConstruction ( IntMat & projBasis,
const IntMat & inBasis,
const Coordinates & proj,
const Int & m,
const ProjConstructType projType = LLLPROJ,
const double delta = 0.9 )

In this version, the construction method is passed as a parameter.

The default is LLL. In the triangular case, a temporary matrix is created internally.

◆ toStringNorm()

static std::string LatticeTester::toStringNorm ( NormType norm)
static

The following are functions for printing the enum constants in this module.

Each function returns the value of the enum variable given as input as a string.

◆ toStringOutput()

static std::string LatticeTester::toStringOutput ( OutputType out)
static

◆ toStringProblem()

static std::string LatticeTester::toStringProblem ( ProblemType prob)
static

◆ toStringPrecision()

static std::string LatticeTester::toStringPrecision ( PrecisionType precision)
static

◆ toStringPrime()

static std::string LatticeTester::toStringPrime ( PrimeType prim)
static

◆ toStringCriterion()

static std::string LatticeTester::toStringCriterion ( CriterionType criter)
static

◆ toStringNorma()

static std::string LatticeTester::toStringNorma ( NormaType norma)
static

◆ toStringCalc()

static std::string LatticeTester::toStringCalc ( CalcType calc)
static

◆ toStringReduction()

static std::string LatticeTester::toStringReduction ( ReductionType reduct)
static

◆ toStringDecomp()

static std::string LatticeTester::toStringDecomp ( DecompTypeBB decomp)
static

◆ toStringProjConstruct()

static std::string LatticeTester::toStringProjConstruct ( ProjConstructType proj)
static

◆ toStringMeritType()

static std::string LatticeTester::toStringMeritType ( MeritType merit)
static

◆ redLLL() [1/7]

template<typename Int , typename Real >
static void LatticeTester::redLLL ( IntMat & basis,
double delta = 0.99999,
long dim = 0,
RealVec * sqlen = 0 )
static

This function uses the NTL implementation of the LLL reduction algorithm with factor delta, presented in [9] (see also [7]).

The reduction is applied to the first dim basis vectors and coordinates when dim > 0, and to the entire basis (all vectors) when dim=0. In the former case, the transformations are not applied to all the columns, so we will no longer have a consistent basis for the original lattice if it had more than dim dimensions. To recover a basis for the full lattice in this case, we may save it before calling this function, or rebuild it.

This function always uses the Euclidean norm. The factor delta must be in [1/2, 1). The closer it is to 1, the more the basis is reduced, in the sense that the LLL algorithm will enforce tighter conditions on the basis. The returned basis always has its shortest vector in first place. The vector pointed by sqlen (if given) will contain the square lengths of the basis vectors. To recover these values in a Vec<double> v one can pass &v to the function. The parameter precision specifies the precision of the floating point numbers that the algorithm will use. EnumTypes.h provides a list of the possible values, and their description is done in the module LLL of NTL.

◆ redLLLExact() [1/3]

template<typename Int >
static void LatticeTester::redLLLExact ( IntMat & basis,
double delta = 0.99999 )
static

This static function implements an exact algorithm from NTL to perform the original LLL reduction.

This is slower than redLLLNTL, but more accurate. It does not take the dim and sqlen parameters (for now).

◆ redBKZ() [1/7]

template<typename Int , typename Real >
static void LatticeTester::redBKZ ( IntMat & basis,
double delta = 0.99999,
int64_t blocksize = 10,
long prune = 0,
long dim = 0,
RealVec * sqlen = 0 )
static

This calls the NTL implementation of the floating point version of the BKZ reduction algorithm presented in [9], with reduction factor delta, block size blocksize, pruning parameter ‘prune’; see [7].

The other paraameters have the same meaning as in redLLLNTL. The parameter blocksize gives the size of the blocks in the BKZ reduction. Roughly, larger blocks means a stronger condition. A blocksize of 2 is equivalent to LLL reduction.

◆ redLLL() [2/7]

template<typename Int , typename Real >
void LatticeTester::redLLL ( IntMat & basis,
double delta = 0.99999,
long dim = 0,
RealVec * sqlen = 0 )

This function uses the NTL implementation of the LLL reduction algorithm with factor delta, presented in [9] (see also [7]).

The reduction is applied to the first dim basis vectors and coordinates when dim > 0, and to the entire basis (all vectors) when dim=0. In the former case, the transformations are not applied to all the columns, so we will no longer have a consistent basis for the original lattice if it had more than dim dimensions. To recover a basis for the full lattice in this case, we may save it before calling this function, or rebuild it.

This function always uses the Euclidean norm. The factor delta must be in [1/2, 1). The closer it is to 1, the more the basis is reduced, in the sense that the LLL algorithm will enforce tighter conditions on the basis. The returned basis always has its shortest vector in first place. The vector pointed by sqlen (if given) will contain the square lengths of the basis vectors. To recover these values in a Vec<double> v one can pass &v to the function. The parameter precision specifies the precision of the floating point numbers that the algorithm will use. EnumTypes.h provides a list of the possible values, and their description is done in the module LLL of NTL.

◆ redLLL() [3/7]

template<>
void LatticeTester::redLLL ( NTL::Mat< int64_t > & basis,
double delta,
long dim,
NTL::Vec< double > * sqlen )

◆ redLLL() [4/7]

template<>
void LatticeTester::redLLL ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long dim,
NTL::Vec< double > * sqlen )

◆ redLLL() [5/7]

template<>
void LatticeTester::redLLL ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long dim,
NTL::Vec< xdouble > * sqlen )

◆ redLLL() [6/7]

template<>
void LatticeTester::redLLL ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long dim,
NTL::Vec< quad_float > * sqlen )

◆ redLLL() [7/7]

template<>
void LatticeTester::redLLL ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long dim,
NTL::Vec< NTL::RR > * sqlen )

◆ redLLLExact() [2/3]

template<typename Int >
void LatticeTester::redLLLExact ( IntMat & basis,
double delta = 0.99999 )

This static function implements an exact algorithm from NTL to perform the original LLL reduction.

This is slower than redLLLNTL, but more accurate. It does not take the dim and sqlen parameters (for now).

◆ redLLLExact() [3/3]

template<>
void LatticeTester::redLLLExact ( NTL::Mat< NTL::ZZ > & basis,
double delta )

◆ redBKZ() [2/7]

template<typename Int , typename Real >
void LatticeTester::redBKZ ( IntMat & basis,
double delta,
long blocksize,
long prune,
long dim,
RealVec * sqlen )

◆ redBKZ() [3/7]

template<>
void LatticeTester::redBKZ ( NTL::Mat< int64_t > & basis,
double delta,
long blocksize,
long prune,
long dim,
NTL::Vec< double > * sqlen )

◆ redBKZ() [4/7]

template<>
void LatticeTester::redBKZ ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long blocksize,
long prune,
long dim,
NTL::Vec< double > * sqlen )

◆ redBKZ() [5/7]

template<>
void LatticeTester::redBKZ ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long blocksize,
long prune,
long dim,
NTL::Vec< xdouble > * sqlen )

◆ redBKZ() [6/7]

template<>
void LatticeTester::redBKZ ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long blocksize,
long prune,
long dim,
NTL::Vec< quad_float > * sqlen )

◆ redBKZ() [7/7]

template<>
void LatticeTester::redBKZ ( NTL::Mat< NTL::ZZ > & basis,
double delta,
long blocksize,
long prune,
long dim,
NTL::Vec< NTL::RR > * sqlen )

◆ swap9()

template<typename T >
void LatticeTester::swap9 ( T & x,
T & y )
inline

Takes references to two variables of a generic type and swaps their content.

This uses the assignment operator, so it might not always work well if this operator's implementation is not thorough.

◆ RandU01()

double LatticeTester::RandU01 ( )

Returns a random number in \([0, 1)\).

The number will have 53 pseudo-random bits.

◆ RandInt()

int64_t LatticeTester::RandInt ( int64_t i,
int64_t j )

Return a uniform pseudo-random integer in \([i, j]\).

Note that the numbers \(i\) and \(j\) are part of the possible output. It is important that \(i < j\) because the underlying arithmetic uses unsigned integers to store j-i+1 and that will be undefined behavior.

◆ RandBits()

std::uint64_t LatticeTester::RandBits ( int64_t s)

Returns the first s pseudo-random bits of the underlying RNG in the form of a s-bit integer.

It is imperative that \(1 \leq s \leq 64\) because the RNG is 64 bits wide.

◆ SetSeed()

void LatticeTester::SetSeed ( std::uint64_t seed)

Sets the seed of the generator.

Because of the constraints on the state, seed has to be \( > 2\). If this is never called, a default seed will be used.

◆ mysqrt()

double LatticeTester::mysqrt ( double x)
inline

Returns \(\sqrt{x}\) for \(x\ge0\), and \(-1\) for \(x < 0\).

◆ Lg() [1/2]

template<typename T >
double LatticeTester::Lg ( const T & x)
inline

Returns the logarithm of \(x\) in base 2.

◆ Lg() [2/2]

double LatticeTester::Lg ( std::int64_t x)
inline

Returns the logarithm of \(x\) in base 2.

◆ abs()

template<typename Scal >
Scal LatticeTester::abs ( Scal x)
inline

Returns the absolute value of \(x\).

◆ sign()

template<typename T >
std::int64_t LatticeTester::sign ( const T & x)
inline

Returns the sign of x.

The sign is 1 if x>0, 0 if x==0 and -1 if x<0.

◆ Round()

template<typename Real >
Real LatticeTester::Round ( Real x)
inline

Returns the value of x rounded to the NEAREST integer value.

(This does not truncate the integer value as is usual in computer arithmetic.)

◆ Factorial()

std::int64_t LatticeTester::Factorial ( int64_t t)

Calculates \(t!\), the factorial of \(t\) and returns it as an std::int64_t.

◆ Quotient()

template<typename Int >
void LatticeTester::Quotient ( const Int & a,
const Int & b,
Int & q )
inline

Computes a/b, truncates the fractional part and puts the result in q.

This function is overloaded to work as specified on NTL::ZZ integers. Example:

\(a\) \(b\) \(q\)
\(5\) 3 1
\(-5\) \(3\) \(-1\)
\(5\) \(-3\) \(-1\)
\(-5\) \(-3\) \(1\)

◆ Modulo()

template<typename Real >
void LatticeTester::Modulo ( const Real & a,
const Real & b,
Real & r )
inline

Computes the remainder of a/b and stores its positive equivalent mod b in r.

This works with std::int64_t, NTL::ZZ and real numbers.

◆ ModuloTowardZero() [1/2]

template<typename Int >
void LatticeTester::ModuloTowardZero ( const Int & a,
const Int & b,
Int & r )
inline

Computes the element r = a mod b that is closest to 0.

Assumes b > 0. For example, 8 mod 11 = -3, -2 mod 11 = -2, -10 mod 11 = 1.

◆ ModuloTowardZero() [2/2]

void LatticeTester::ModuloTowardZero ( const std::int64_t & a,
const std::int64_t & b,
std::int64_t & r )
inline

◆ Divide()

template<typename Real >
void LatticeTester::Divide ( Real & q,
Real & r,
const Real & a,
const Real & b )
inline

Computes the quotient \(q = a/b\) and remainder \(r = a \bmod b\).

Truncates \(q\) to the nearest integer toward 0. One always has \(a = qb + r\) and \(|r| < |b|\). This works with std::int64_t, NTL::ZZ and real numbers.

\(a\) \(b\) \(q\) \(r\)
\(5\) 3 1 \(2\)
\(-5\) \(3\) \(-1\) \(-2\)
\(5\) \(-3\) \(-1\) \(2\)
\(-5\) \(-3\) \(1\) \(-2\)

◆ DivideRound()

template<typename Real >
void LatticeTester::DivideRound ( const Real & a,
const Real & b,
Real & q )
inline

Computes \(a/b\), rounds the result to the nearest integer and returns the result in \(q\).

This works with std::int64_t, NTL::ZZ and real numbers.

◆ gcd()

std::int64_t LatticeTester::gcd ( std::int64_t a,
std::int64_t b )

Returns the value of the greatest common divisor of \(a\) and \(b\) by using Stein's binary GCD algorithm.

◆ Euclide() [1/2]

template<typename Int >
void LatticeTester::Euclide ( const Int & A,
const Int & B,
Int & C,
Int & D,
Int & E,
Int & F,
Int & G )

This method computes the greater common divisor of A and B with Euclid's algorithm.

This will store this gcd in G and also the linear combination that permits to get G from A and B. This function should work with std::int64_t and NTL::ZZ.

For \(A\) and \(B\) this will assign to \(C\), \(D\), \(E\), \(F\) and \(G\) values such that:

\begin{align*} C a + D b & = G = \mbox{GCD } (a,b)\\ E a + F b & = 0. \end{align*}

◆ Euclide() [2/2]

template<typename Int >
void LatticeTester::Euclide ( const Int & A,
const Int & B,
Int & C,
Int & D,
Int & G )

◆ TransposeMatrix()

template<typename Int >
void LatticeTester::TransposeMatrix ( IntMat & mat,
IntMat & mat2 )

Takes a set of generating vectors in the matrix mat and iteratively transforms it into an upper triangular lattice basis into the matrix mat2.

mat and mat2 have to have the same number of rows and the same number of columns. All the computations will be done modulo mod, which means that you must know the rescaling factor for the vector system to call this function. After the execution, mat will be a matrix containing irrelevant information and mat2 will contain an upper triangular basis.

For more details please look at [latTesterGide]. This algorithm basically implements what is written in this guide. The matrix mat contains the set of vectors that is used and modified at each step to get a new vector from the basis. Lower triangularization

◆ CreateVect() [1/2]

template<typename Real >
void LatticeTester::CreateVect ( Real *& A,
int64_t d )
inline

Allocates memory to A as an array of Real of dimension d and initializes its elements to 0.

Real has to be a numeric type.

◆ CreateVect() [2/2]

template<typename Vect >
void LatticeTester::CreateVect ( Vect & A,
int64_t d )
inline

Creates the vector A of dimensions d+1 and initializes its elements to 0.

The type Vect must have a method SetLength, as for Vec<T> in NTL.

◆ DeleteVect() [1/2]

template<typename Real >
void LatticeTester::DeleteVect ( Real *& A)
inline

Frees the memory used by the vector A.

This calls delete[] on A so trying to access A after using this is unsafe.

◆ DeleteVect() [2/2]

template<typename Vect >
void LatticeTester::DeleteVect ( Vect & A)
inline

Frees the memory used by the vector A, destroying all the elements it contains.

Vect type has to have a kill() method that deallocates all the elements in the vector.

◆ SetZero() [1/2]

template<typename Real >
void LatticeTester::SetZero ( Real * A,
int64_t d )
inline

Sets the first d of A to 0.

◆ SetZero() [2/2]

template<typename Vect >
void LatticeTester::SetZero ( Vect & A,
int64_t d )
inline

Sets the first d components of A to 0.

◆ SetValue()

template<typename Real >
void LatticeTester::SetValue ( Real * A,
int64_t d,
const Real & x )
inline

Sets the first d components of A to 0.

Sets the first d components of A to the value x.

◆ toString() [1/2]

template<typename Vect >
std::string LatticeTester::toString ( const Vect & A,
int64_t c,
int64_t d,
const char * sep = " " )

Returns a string containing A[c] to A[d-1] formated as [A[c]sep...sepA[d-1]].

In this string, components are separated by string sep. By default, sep is just a whitespace character.

◆ toString() [2/2]

template<typename Vect >
std::string LatticeTester::toString ( const Vect & A,
int64_t d )

Returns a string containing the first d components of the vector A as a string.

Calls toString(const Vect&, int, int, const char*).

◆ ProdScal()

template<typename Int , typename Real >
void LatticeTester::ProdScal ( const IntVec & A,
const IntVec & B,
int64_t n,
Real & D )
inline

Computes the scalar product of vectors A and B truncated to their n first components, then puts the result in D.

There is a lot to consider when passing types to this function. The best is for Vect1 to be the same type as Vect2 and for Scal to be the same as Int, and that those types are the ones stored in Vect1 and Vect2.

WARNING: This uses so many types without check about them and also assumes all those types can be converted to each other without problem. This is used in some places to compute a floating point norm of vectors with integers values. Take care when using this function.

◆ Invert()

template<typename Int >
void LatticeTester::Invert ( const IntVec & A,
IntVec & B,
int64_t n )
inline

Takes an input vector A of dimension n+1 and fill the vector B with the values [-A[n] -A[n-1] ... -A[1][1].

B is assumed to be of dimension at least n+1.

◆ CalcNorm()

template<typename Int , typename Real >
void LatticeTester::CalcNorm ( const IntVec & V,
int64_t n,
Real & S,
NormType norm )
inline

Computes the norm norm of vector V trunctated to its n first components, and puts the result in S.

Scal has to be a floating point type. For the L2 norm, it returns the square norm instead.

◆ CopyPartVec()

template<typename Vect >
void LatticeTester::CopyPartVec ( Vect & toVec,
const Vect & fromVec,
int64_t c )
inline

Copies the first c components of vector fromVec into vector toVec.

◆ CopyPartMat()

template<typename Matr >
void LatticeTester::CopyPartMat ( Matr & toMat,
const Matr & fromMat,
int64_t r,
int64_t c )
inline

Copies the first r rows and c columns of matrix fromMat into matrix toMat.

◆ CopyVect()

template<typename Vect >
void LatticeTester::CopyVect ( Vect & A,
const Vect & B,
int64_t n )
inline

Copies the first n components of vector B into vector A.

◆ ModifVect()

template<typename Vect , typename Scal >
void LatticeTester::ModifVect ( Vect & A,
const Vect & B,
Scal x,
int64_t n )
inline

Adds the first n components of vector B multiplied by x to first n components of vector A.

This will modify A. This does weird type conversions and may not work well if different types are used.

◆ ModifVectModulo()

template<typename Vect , typename Scal , typename Int >
void LatticeTester::ModifVectModulo ( Vect & A,
const Vect & B,
Scal x,
Int m,
int64_t n )
inline

Adds the first n components of vector B multiplied by x to first n components of vector A, modulo m.

This will modify A. The elements that are not multiples of m are reduced mod m.

◆ ChangeSign()

template<typename Vect >
void LatticeTester::ChangeSign ( Vect & A,
int64_t n )
inline

Changes the sign (multiplies by -1) the first n components of vector A.

◆ GCD2vect()

std::int64_t LatticeTester::GCD2vect ( std::vector< std::int64_t > V,
int64_t k,
int64_t n )
inline

Computes the greatest common divisor of V[k],...,V[n-1].

◆ CreateMatr() [1/4]

template<typename Real >
void LatticeTester::CreateMatr ( Real **& A,
int64_t d )
inline

Allocates memory to a square matrix A of dimensions \(d \times d\) and initializes its elements to 0.

◆ CreateMatr() [2/4]

template<typename Real >
void LatticeTester::CreateMatr ( Real **& A,
int64_t line,
int64_t col )
inline

Allocates memory for the matrix A of dimensions \(\text{line} \times \text{col}\) and initializes its elements to 0.

◆ CreateMatr() [3/4]

template<typename Int >
void LatticeTester::CreateMatr ( IntMat & A,
int64_t d )
inline

Resizes the matrix A to a square matrix of dimensions d*d and re-initializes its elements to 0.

◆ CreateMatr() [4/4]

template<typename Int >
void LatticeTester::CreateMatr ( IntMat & A,
int64_t line,
int64_t col )
inline

Resizes the matrix A to a matrix of dimensions \(line \times col\) and re-initializes its elements to 0.

◆ DeleteMatr() [1/3]

template<typename Real >
void LatticeTester::DeleteMatr ( Real **& A,
int64_t d )
inline

Frees the memory used by the \(d \times d\) matrix A.

This will not free all the memory allocated to A if A is of greater dimension and it can cause a memory leak.

◆ DeleteMatr() [2/3]

template<typename Real >
void LatticeTester::DeleteMatr ( Real **& A,
int64_t line,
int64_t col )
inline

Frees the memory used by the matrix A of dimension \(\text{line} \times \text{col}\).

This will not free all the memory allocated to A if A is of greater dimension and it can cause a memory leak.

◆ DeleteMatr() [3/3]

template<typename Int >
void LatticeTester::DeleteMatr ( IntMat & A)
inline

Calls the clear() method on A.

A has to have a clear() method that frees the memory allocated to it.

◆ CopyMatr() [1/2]

template<typename Matr >
void LatticeTester::CopyMatr ( Matr & A,
const Matr & B,
int64_t n )
inline

Copies the \(n \times n\) submatrix of the first lines and columns of B into matrix A.

This function does not check for sizes, so A and B both have to be at leat \(n \times n\).

◆ CopyMatr() [2/2]

template<typename Matr >
void LatticeTester::CopyMatr ( Matr & A,
const Matr & B,
int64_t line,
int64_t col )
inline

Copies the \(\text{line} \times col\) submatrix of the first lines and columns of B into matrix A.

This function does not check for sizes, so A and B both have to be at leat \(line \times col\).

◆ toStr()

template<typename MatT >
std::string LatticeTester::toStr ( const MatT & mat,
int64_t d1,
int64_t d2,
int64_t prec = 2 )

Returns a string that is a representation of mat.

This string represents the \(d1 \times d2\) submatrix of the first lines and columns of mat.

◆ ProductDiagonal()

template<typename Int >
void LatticeTester::ProductDiagonal ( const NTL::Mat< Int > & A,
long dim,
Int & prod )

Returns the product of the diagonal elements of the matrix A, which is assumed to be square dim x dim.

◆ CheckTriangular()

template<typename Int >
bool LatticeTester::CheckTriangular ( const NTL::Mat< Int > & A,
long dim,
const Int m )

Checks that the upper \(\text{dim} \times \text{dim}\) submatrix of A is triangular modulo m.

This will return true if all the elements under the diagonal are equal to zero modulo m and false otherwise. If m is 0, this function simply verifies that the matrix is triangular.

◆ checkInverseModm()

template<typename Int >
bool LatticeTester::checkInverseModm ( const NTL::Mat< Int > & A,
const NTL::Mat< Int > & B,
const Int m )

Checks if A and B are m-dual to each other.

They must be square with the same dimensions. Returns true if AB = mI, false otherwise.

◆ calcDual()

template<typename Matr , typename Int >
void LatticeTester::calcDual ( const Matr & A,
Matr & B,
int64_t d,
const Int & m )

Takes an upper triangular basis A and computes an m-dual lattice basis to this matrix.

For this algorithm to work, A has to be upper triangular and all the coefficients on the diagonal have to divide m.

For B to be m-dual to A, we have to have that \(AB^t = mI\). It is quite easy to show that, knowing A is upper triangular, B will be a lower triangular matrix with A(i,i)*B(i,i) = m for all i and \( A_i \cdot B_j = 0\) for \(i\neq j\). To get the second condition, we simply have to recursively take for each line

\[B_{i,j} = -\frac{1}{A_{j,j}}\sum_{k=j+1}^i A_{j,k} B_{i,k}.\]

◆ myExit()

void LatticeTester::myExit ( std::string msg)
inline

Simple error exit function, prints msg on exit.

◆ operator<<() [1/6]

template<class K , class T , class C , class A >
std::ostream & LatticeTester::operator<< ( std::ostream & out,
const std::map< K, T, C, A > & x )

Streaming operator for maps.

Formats a map as: { key1=>val1, ..., keyN=>valN }.

◆ operator<<() [2/6]

template<class T1 , class T2 >
std::ostream & LatticeTester::operator<< ( std::ostream & out,
const std::pair< T1, T2 > & x )

Streaming operator for vectors.

Formats a pair as: (first,second).

◆ operator<<() [3/6]

template<class T , class A >
std::ostream & LatticeTester::operator<< ( std::ostream & out,
const std::vector< T, A > & x )

Streaming operator for vectors.

Formats a vector as: [ val1, ..., valN ].

◆ operator<<() [4/6]

template<class K , class C , class A >
std::ostream & LatticeTester::operator<< ( std::ostream & out,
const std::set< K, C, A > & x )

Streaming operator for sets.

Formats a set as: { val1, ..., valN }.

◆ copyMatrixToMat()

template<typename Matr1 , typename Matr2 >
void LatticeTester::copyMatrixToMat ( Matr1 & A,
Matr2 & B )

◆ printBase()

template<typename Int >
void LatticeTester::printBase ( IntMat bas_mat)

◆ printBase2()

template<typename Int >
void LatticeTester::printBase2 ( IntMat bas_mat)

◆ copy() [1/2]

template<typename Int >
void LatticeTester::copy ( IntMat & b1,
IntMat & b2 )

◆ copy() [2/2]

template<typename Int >
void LatticeTester::copy ( IntMat & b1,
IntMat & b2,
int64_t r,
int64_t c )

◆ getWidth()

int64_t LatticeTester::getWidth ( clock_t time[],
int64_t dim,
std::string message,
clock_t totals[],
int64_t ind )
inline

Variable Documentation

◆ MAX_LONG_DOUBLE

const double LatticeTester::MAX_LONG_DOUBLE = 9007199254740992.0

Maximum integer that can be represented exactly as a double: \(2^{53}\).

◆ TWO_EXP

const std::int64_t LatticeTester::TWO_EXP[]
extern

Table of powers of 2: TWO_EXP[ \(i\)] \(= 2^i\), \(i=0, 1, …, 63\).