Lattice Tester Online Documentation 0.1.0-861
Software Package For Testing The Uniformity Of Integral Lattices In The Real Space
|
Lattice namespace. More...
Namespaces | |
namespace | CoordinateSets |
The classes FromRange , SubSets , and AddCoordinate are defined here. | |
namespace | Random |
This class generates random numbers (in fact pseudo-random numbers). | |
Classes | |
class | Chrono |
This class acts as an interface to the system clock to compute the CPU time used by parts of a program. More... | |
class | Coordinates |
This is basically a std::set<std::size_t> , i.e., a standard C++ set. More... | |
class | FigureOfMeritDualM |
This class offers funcitons to calculate the figure of merit (FOM) for the m-dual of any given 'IntLatticeExt' object. More... | |
class | FigureOfMeritM |
$\f, for each order \(s > 1\). More... | |
class | FoMCalc |
This class offers methods (functions) to calculate figure of merit for a given IntLattice 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 | Lacunary |
REMOVE! More... | |
class | MRGLattice |
This subclass of IntLatticeExt defines an MRG lattice and is similar to Rank1Lattice, but constructs lattices associated with multiple recursive generators (MRGs) with modulus m, order k, and vector of multipliers a = (a_1, . More... | |
class | NormaBestBound |
In this normalizer, the Hermite constants \(\gamma_s\) are approximated using the best upper bounds that are available. More... | |
class | NormaBestLat |
This Normalizer class implements upper bounds on the length of the shortest nonzero vector in a lattice. More... | |
class | NormaLaminated |
This Normalizer class implements 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 | NormaMinkL2 |
This class implements Minkowski’s theoretical LOWER bound on the length of the shortest non-zero vector in a lattice, with the \({\mathcal{L}}_2\) norm. More... | |
class | NormaPalpha |
This class implements some 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 | ParamReader |
Utility class that can be used to read different kind of data from a file. More... | |
class | Rank1Lattice |
This subclass of IntLatticeExt defines a general rank 1 lattice rule in \(d\) dimensions, whose points \(\mathbb{u}_i\) are defined by. More... | |
class | ReducerBB |
This ReducerBB class provides facilities to reduce the basis of a lattice (an IntLattice object) in various ways (pairwise, LLL, BKZ, Minkowski [4], [9], [10]), and to find a shortest nonzero vector in the lattice using a BB algorithm [rFIN85a]. More... | |
class | Types |
Sets standard typedef ’s for the types that can be used in LatticeTester. 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. More... | |
class | WeightsPOD |
Defines product and order-dependent (POD) weights. More... | |
class | WeightsProduct |
Defines product weights. More... | |
class | WeightsProjectionDependent |
Defines projection-dependent weights. More... | |
class | WeightsUniform |
Specifies projection weights that are the same (usually 1) for all projections. More... | |
class | Writer |
This is an abstract class that represents an interface to Writer classes. More... | |
class | WriterRes |
This class is a simple implementation of the Writer abstract class to write in plain text format on the stream. 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 , BESTBOUND , LAMINATED , ROGERS , MINKL1 , MINKL2 , 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 IntMat , typename RealVec > | |
static long | LLLConstruction0 (IntMat &gen, const double delta=0.9, long r=0, long c=0, RealVec *sqlen=0) |
This static class offers methods (functions) to construct a basis from a set of generating vectors that are not necessarily independent, to construct a triangular basis, to construct the basis for a projection over a given subset of coordinates, and to obtain the \(m\)-dual of a given basis. | |
template<typename IntMat , typename Int , typename RealVec > | |
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 in case the set of generating vectors do not generate a full-dimensional lattice, it adds the vectors \(m e_i\) to the generating set, so it always returns a square matrix. | |
template<typename IntMat , typename Int > | |
static void | lowerTriangularBasis (IntMat &gen, IntMat &basis, 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 IntMat , typename Int > | |
static void | upperTriangularBasis (IntMat &gen, IntMat &basis, const Int &m, long r=0, long c=0) |
Same as lowerTriangularBasis , except that the returned basis is upper triangular. | |
template<typename IntMat , typename Int > | |
static void | mDualUpperTriangular (const IntMat &basis, IntMat &basisDual, const Int &m, long dim=0) |
Takes an upper triangular basis matrix basis and computes the m-dual basis basisDual . | |
template<typename IntMat , typename Int > | |
static void | mDualUpperTriangular96 (IntMat &basis, IntMat &basisDual, const Int &m, long dim=0) |
This function does essentially the same thing as mDualUpperTriangular , but the algorithm is slightly different. | |
template<typename IntMat , typename Int > | |
static void | mDualBasis (const IntMat &basis, IntMat &basisDual, 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 IntMat > | |
static void | projectMatrix (const IntMat &in, IntMat &out, 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 IntMat , typename Int , typename RealVec > | |
static void | projectionConstructionLLL (const IntMat &inBasis, IntMat &projBasis, 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 IntMat , typename Int > | |
static void | projectionConstructionUpperTri (const IntMat &inBasis, IntMat &projBasis, 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 IntMat , typename Int > | |
static void | projectionConstructionUpperTri (const IntMat &inBasis, IntMat &projBasis, const Coordinates &proj, const Int &m, long r=0) |
template<typename IntMat , typename Int > | |
static void | projectionConstruction (const IntMat &inBasis, IntMat &projBasis, 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::matrix< long > &gen, const double delta, long r, long c, NTL::vector< double > *sqlen) |
template<> | |
long | LLLConstruction0 (NTL::matrix< NTL::ZZ > &gen, const double delta, long r, long c, NTL::vector< double > *sqlen) |
template<> | |
long | LLLConstruction0 (NTL::matrix< NTL::ZZ > &gen, const double delta, long r, long c, NTL::vector< xdouble > *sqlen) |
template<> | |
long | LLLConstruction0 (NTL::matrix< NTL::ZZ > &gen, const double delta, long r, long c, NTL::vector< quad_float > *sqlen) |
template<> | |
long | LLLConstruction0 (NTL::matrix< NTL::ZZ > &gen, const double delta, long r, long c, NTL::vector< NTL::RR > *sqlen) |
template<typename IntMat , typename Int , typename RealVec > | |
void | LLLBasisConstruction (IntMat &gen, const Int &m, double delta, long r, long c, RealVec *sqlen) |
Similar to LLLConstruction0 , except that in case the set of generating vectors do not generate a full-dimensional lattice, it adds the vectors \(m e_i\) to the generating set, so it always returns a square matrix. | |
template<typename IntMat , typename Int > | |
void | lowerTriangularBasis (IntMat &gen, IntMat &basis, 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 IntMat , typename Int > | |
void | upperTriangularBasis (IntMat &gen, IntMat &basis, const Int &m, long dim1, long dim2) |
Same as lowerTriangularBasis , except that the returned basis is upper triangular. | |
template<typename IntMat , typename Int > | |
void | mDualUpperTriangular (const IntMat &A, IntMat &B, const Int &m, long dim) |
For B to be m -dual to A , we have to have that \(AB^t = mI\). | |
template<typename IntMat , typename Int > | |
void | mDualUpperTriangular96 (IntMat &basis, IntMat &basisDual, const Int &m, long dim) |
This function does essentially the same thing as mDualUpperTriangular , but the algorithm is slightly different. | |
template<> | |
void | mDualUpperTriangular96 (NTL::matrix< NTL::ZZ > &basis, NTL::matrix< NTL::ZZ > &basisDual, const NTL::ZZ &m, long dim) |
template<typename IntMat , typename Int > | |
void | mDualBasis (const IntMat &basis, IntMat &basisDual, 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<> | |
void | mDualBasis (const NTL::matrix< NTL::ZZ > &basis, NTL::matrix< NTL::ZZ > &basisDual, const NTL::ZZ &m) |
template<typename IntMat > | |
void | projectMatrix (const IntMat &in, IntMat &out, 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 IntMat , typename Int , typename RealVec > | |
void | projectionConstructionLLL (const IntMat &inBasis, IntMat &projBasis, 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 IntMat , typename Int > | |
void | projectionConstructionUpperTri (const IntMat &inBasis, IntMat &projBasis, 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 IntMat , typename Int > | |
void | projectionConstructionUpperTri (const IntMat &inBasis, IntMat &projBasis, const Coordinates &proj, const Int &m, long r) |
template<typename IntMat , typename Int > | |
void | projectionConstruction (const IntMat &inBasis, IntMat &projBasis, const Coordinates &proj, const Int &m, const ProjConstructType projType, const double delta) |
In this version, the construction method is passed as a parameter. | |
std::string | toString (Chrono &timer) |
Returns the value of the duration from timer . | |
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) |
double | sqrtReal (const double &a) |
xdouble | sqrtReal (const xdouble &a) |
quad_float | sqrtReal (const quad_float &a) |
NTL::RR | sqrtReal (const NTL::RR &a) |
std::int64_t | lFactorial (int64_t t) |
Calculates \(t!\), the factorial of \(t\). | |
double | Digamma (double x) |
Returns the value of the logarithmic derivative of the Gamma function \(\psi(x) = \Gamma'(x) / \Gamma(x)\). | |
double | BernoulliPoly (int64_t n, double x) |
Evaluates the Bernoulli polynomial \(B_n(x)\) of degree \(n\) at \(x\). | |
double | Harmonic (std::int64_t n) |
Computes the \(n\)-th harmonic number \(H_n = \sum_{j=1}^n 1/j\). | |
double | Harmonic2 (std::int64_t n) |
Computes the sum. | |
double | FourierC1 (double x, std::int64_t n) |
Computes and returns the value of the series (see [5]) | |
double | FourierE1 (double x, std::int64_t n) |
Computes and returns the value of the series. | |
void | negativeCholesky () |
template<typename IntMat , typename RealVec > | |
static void | redLLL (IntMat &basis, double delta=0.99999, long dim=0, RealVec *sqlen=0) |
This file provides only static functions to reduce a lattice basis in some ways (pairwise, LLL, BKZ [4], [9], [10]). | |
template<typename IntMat > | |
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 IntMat , typename RealVec > | |
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 [10], with reduction factor delta , block size blocksize , pruning parameter ‘prune’; see [7]. | |
template<typename IntMat , typename RealVec > | |
void | redLLL (IntMat &basis, double delta, long dim, RealVec *sqlen) |
This file provides only static functions to reduce a lattice basis in some ways (pairwise, LLL, BKZ [4], [9], [10]). | |
template<> | |
void | redLLL (NTL::matrix< int64_t > &basis, double delta, long dim, NTL::vector< double > *sqlen) |
template<> | |
void | redLLL (NTL::matrix< NTL::ZZ > &basis, double delta, long dim, NTL::vector< double > *sqlen) |
template<> | |
void | redLLL (NTL::matrix< NTL::ZZ > &basis, double delta, long dim, NTL::vector< xdouble > *sqlen) |
template<> | |
void | redLLL (NTL::matrix< NTL::ZZ > &basis, double delta, long dim, NTL::vector< quad_float > *sqlen) |
template<> | |
void | redLLL (NTL::matrix< NTL::ZZ > &basis, double delta, long dim, NTL::vector< NTL::RR > *sqlen) |
template<typename IntMat > | |
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::matrix< NTL::ZZ > &basis, double delta) |
template<typename IntMat , typename RealVec > | |
void | redBKZ (IntMat &basis, double delta, long blocksize, long prune, long dim, RealVec *sqlen) |
template<> | |
void | redBKZ (NTL::matrix< int64_t > &basis, double delta, long blocksize, long prune, long dim, NTL::vector< double > *sqlen) |
template<> | |
void | redBKZ (NTL::matrix< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::vector< double > *sqlen) |
template<> | |
void | redBKZ (NTL::matrix< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::vector< xdouble > *sqlen) |
template<> | |
void | redBKZ (NTL::matrix< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::vector< quad_float > *sqlen) |
template<> | |
void | redBKZ (NTL::matrix< NTL::ZZ > &basis, double delta, long blocksize, long prune, long dim, NTL::vector< 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 IntMat > | |
void | printBase (IntMat bas_mat) |
template<typename IntMat > | |
void | printBase2 (IntMat bas_mat) |
template<typename IntMat > | |
void | copy (IntMat &b1, IntMat &b2) |
template<typename IntMat > | |
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) |
double | Harmonic (long n) |
double | Harmonic2 (long n) |
double | FourierC1 (double x, long n) |
double | FourierE1 (double x, long n) |
std::int64_t | RandInt (std::int64_t i, std::int64_t j) |
static bool | check_next_chars (istream &is, const string &token) |
Helper function to check the next characters from an input stream. | |
static void | skip_any (istream &is, const string &characters) |
Helper function to skip all characters of a given class. | |
istream & | operator>> (std::istream &is, WeightsProjectionDependent &weights) |
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 | |
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 fractionnal 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 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 IntMat > | |
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. | |
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 Vect1 , typename Vect2 , typename Scal > | |
void | ProdScal (const Vect1 &A, const Vect2 &B, int64_t n, Scal &D) |
Computes the scalar product of vectors A and B truncated to their n first components, then puts the result in D . | |
template<typename IntVec > | |
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 Vect , typename Scal > | |
void | CalcNorm (const Vect &V, int64_t n, Scal &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 Vect1 , typename Vect2 , typename Scal > | |
void | ModifVect (Vect1 &A, const Vect2 &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 > | |
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 IntMat > | |
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 IntMat > | |
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 IntMat > | |
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::matrix< 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::matrix< 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::matrix< Int > &A, const NTL::matrix< Int > &B, const Int m) |
Checks if A and B are m-dual to each other. | |
template<typename Matr , typename Int > | |
void | Triangularization (Matr &W, Matr &V, int64_t lin, int64_t col, const Int &m) |
Takes a set of generating vectors in the matrix W and iteratively transforms it into an upper triangular lattice basis into the matrix V . | |
template<typename Matr , typename Int > | |
void | calcDual (const Matr &A, Matr &B, int64_t d, const Int &m) |
Takes a basis A and computes an m-dual lattice basis B. | |
Debugging functions | |
void | MyExit (int64_t status, std::string msg) |
Special exit function. | |
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\). | |
Lattice namespace.
typedef double LatticeTester::Weight |
A scalar weight type.
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 |
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 [11].
???
Enumerator | |
---|---|
LENGTH | |
SPECTRAL | |
BEYER | |
PALPHA | |
BOUND_JS |
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.
BESTBOUND
: 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.
MINKL2
: 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 | |
BESTBOUND | |
LAMINATED | |
ROGERS | |
MINKL1 | |
MINKL2 | |
NONE |
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 |
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 |
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 |
Types of problems that LatticeTester can handle.
Not sure if we still need this. *******
Enumerator | |
---|---|
BASIS | |
DUAL | |
REDUCTION | |
SHORTEST | |
MERIT |
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 |
|
inline |
Returns the absolute value of \(x\).
double LatticeTester::BernoulliPoly | ( | int64_t | n, |
double | x ) |
Evaluates the Bernoulli polynomial \(B_n(x)\) of degree \(n\) at \(x\).
The first Bernoulli polynomials are:
\begin{align*} B_0(x) &= 1 \\ B_1(x) &= x - 1/2 \\ B_2(x) &= x^2-x+1/6 \\ B_3(x) &= x^3 - 3x^2/2 + x/2 \\ B_4(x) &= x^4-2x^3+x^2-1/30 \\ B_5(x) &= x^5 - 5x^4/2 + 5x^3/3 - x/6 \\ B_6(x) &= x^6-3x^5+5x^4/2-x^2/2+1/42 \\ B_7(x) &= x^7 - 7x^6/2 + 7x^5/2 - 7x^3/6 + x/6 \\ B_8(x) &= x^8-4x^7+14x^6/3 - 7x^4/3 +2x^2/3-1/30. \end{align*}
Only degrees \(n \le 8\) are programmed for now.
void LatticeTester::calcDual | ( | const Matr & | A, |
Matr & | B, | ||
int64_t | d, | ||
const Int & | m ) |
Takes a basis A
and computes an m-dual lattice basis B.
The matrix B is the m-dual basis of A. 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}.\]
|
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.
|
inline |
Changes the sign (multiplies by -1) the first n
components of vector A
.
|
static |
Helper function to check the next characters from an input stream.
Returns true
if the next characters in is
are token
. Upon a match, the characters are removed from is
; otherwise, they are left in is
.
bool LatticeTester::checkInverseModm | ( | const NTL::matrix< Int > & | A, |
const NTL::matrix< 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.
bool LatticeTester::CheckTriangular | ( | const NTL::matrix< 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.
void LatticeTester::copy | ( | IntMat & | b1, |
IntMat & | b2 ) |
void LatticeTester::copy | ( | IntMat & | b1, |
IntMat & | b2, | ||
int64_t | r, | ||
int64_t | c ) |
|
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\).
|
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\).
void LatticeTester::copyMatrixToMat | ( | Matr1 & | A, |
Matr2 & | B ) |
|
inline |
Copies the first r
rows and c
columns of matrix fromMat
into matrix toMat
.
|
inline |
Copies the first c
components of vector fromVec
into vector toVec
.
|
inline |
Copies the first n
components of vector B
into vector A
.
|
inline |
Resizes the matrix A
to a square matrix of dimensions d*d
and re-initializes its elements to 0.
|
inline |
Resizes the matrix A
to a matrix of dimensions \(line \times col\) and re-initializes its elements to 0.
|
inline |
Allocates memory to a square matrix A
of dimensions \(d \times d\) and initializes its elements to 0.
|
inline |
Allocates memory for the matrix A
of dimensions \(\text{line} \times
\text{col}\) and initializes its elements to 0.
|
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.
|
inline |
Creates the vector A
of dimensions d+1
and initializes its elements to 0.
The type Vect
has to have a resize(integer_type)
method that sets the size of the instance to the value of the argument.
|
inline |
Calls the clear()
method on A
.
A
has to have a clear()
method that frees the memory allocated to it.
|
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.
|
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.
|
inline |
Frees the memory used by the vector A
.
This calls delete[]
on A
so trying to access A
after using this is unsafe.
|
inline |
Frees the memory used by the vector A
, destroying all the elements it contains.
Vect
type has to have a clear()
method that deallocates all the elements in the vector.
double LatticeTester::Digamma | ( | double | x | ) |
Returns the value of the logarithmic derivative of the Gamma function \(\psi(x) = \Gamma'(x) / \Gamma(x)\).
|
inline |
Computes the quotient \(q = a/b\) and remainder \(r = a \bmod b\).
Truncates \(q\) to the nearest integer towards 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\) |
|
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.
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*}
void LatticeTester::Euclide | ( | const Int & | A, |
const Int & | B, | ||
Int & | C, | ||
Int & | D, | ||
Int & | G ) |
std::int64_t LatticeTester::Factorial | ( | int64_t | t | ) |
Calculates \(t!\), the factorial of \(t\) and returns it as an std::int64_t.
double LatticeTester::FourierC1 | ( | double | x, |
long | n ) |
double LatticeTester::FourierC1 | ( | double | x, |
std::int64_t | n ) |
Computes and returns the value of the series (see [5])
\[ S(x, n) = \sum_{j=1}^{n} \frac{\cos(2\pi j x)}{j}. \]
Restrictions: \(n>0\) and \(0 \le x \le 1\).
double LatticeTester::FourierE1 | ( | double | x, |
long | n ) |
double LatticeTester::FourierE1 | ( | double | x, |
std::int64_t | n ) |
Computes and returns the value of the series.
\[ G(x, n) = \sideset{}{'}\sum_{-n/2<h\le n/2}\; \frac{e^{2\pi i h x}}{|h|}, \]
where the symbol \(\sum^\prime\) means that the term with \(h=0\) is excluded from the sum, and assuming that the imaginary part of \(G(x, n)\) vanishes. Restrictions: \(n>0\) and \(0 \le x \le 1\).
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.
|
inline |
Computes the greatest common divisor of V[k],...,V[n-1]
.
|
inline |
double LatticeTester::Harmonic | ( | long | n | ) |
double LatticeTester::Harmonic | ( | std::int64_t | n | ) |
Computes the \(n\)-th harmonic number \(H_n = \sum_{j=1}^n 1/j\).
double LatticeTester::Harmonic2 | ( | long | n | ) |
double LatticeTester::Harmonic2 | ( | std::int64_t | n | ) |
Computes the sum.
\[ \sideset{}{'}\sum_{-n/2<j\le n/2}\; \frac 1{|j|}, \]
where the symbol \(\sum^\prime\) means that the term with \(j=0\) is excluded from the sum.
|
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
.
std::int64_t LatticeTester::lFactorial | ( | int64_t | t | ) |
Calculates \(t!\), the factorial of \(t\).
Might throw if t
is too large or if std::int64_t can't contain the factorial asked for.
|
inline |
Returns the logarithm of \(x\) in base 2.
|
inline |
Returns the logarithm of \(x\) in base 2.
|
static |
Similar to LLLConstruction0
, except that in case the set of generating vectors do not generate a full-dimensional lattice, it adds the vectors \(m 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.
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 in case the set of generating vectors do not generate a full-dimensional lattice, it adds the vectors \(m 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.
|
static |
This static class offers methods (functions) to construct a basis from a set of generating vectors that are not necessarily independent, to construct a triangular basis, to construct the basis for a projection over a given subset of coordinates, and to obtain the \(m\)-dual of a given basis.
The implementation relies on NTL and uses NTL matrices. When the basis turns out to have fewer rows than columns, some of the methods add implicitly the rescaled unit vectors to the set of generating vectors. In that case, the basis matrix is always square and all the vectors of the form \(m \be_i\) belong to the lattice.
NTL already offers an efficient method to construct an LLL-reduced basis from a set of generating vectors. This is the most effective way of constructing a basis and it is encapsulated in the LLLConstruction0
method given below. This method does not assume that the rescaled unit vectors \(m \be_i\) belong to the lattices and it does not even know about \(m\). The method LLLBasisConstruction
adds those vectors to the set of generating vectors, so it always returns a square basis.
We also offer an alternative methods that construct a triangular basis from a set of generating vectors. They always add the rescaled unit vectors implicitly to the set. The method lowerTriangularBasis
constructs a lower-triangular basis, while upperTriangularBasis
constructs an upper-triangular basis.
To compute the \(m\)-dual of a given basis, we have a general (but slow) method implemented in mDualBasis
, and a much faster method in mDualUpperTriangular
that works only when the basis is upper-triangular.
We also have functions to compute the basis of a projection of a given lattice over a specified set of coordinates. The function projectionConstructionLLL
does this by using LLL to construct the basis of the projection, while projectionConstructionUpperTri
constructs an upper-triangular basis for the projection. The function projectionConstruction
takes the construction method as a parameter.
All functions take as input a IntMat
object that contains either a basis or a set of generating vectors. By default, all rows and columns of this matrix object are used. But in most functions, the user can ask to use only a subset of these rows and columns, via the optional parameters r
and c
. This can permit one to use the same IntMat
object for several numbers of dimensions, to avoid doing many object creations or resizing.
All functions in this class are static, so there is no reason to create any BasisConstruction
object. We also avoid to create new objects (such as vectors and matrices) inside these functions. These functions can be called thousands or millions of times in a program, and we want the user to be able to re-use the same vectors and matrices over and over again instead of creating new ones.
Note that when one of these functions is used for an IntMat
object that store the primal or \(m\)-dual basis inside an IntLattice
object, only the primal or the \(m\)-dual basis is changed, and therefore the \(m\)-dual relationship will no longer hold after the call. The user must be aware of that. Most of the time, there is no need to reinstate this relationship, because for example if we want to compute a shortest nonzero vector is the \(m\)-dual lattice, then once we have an \(m\)-dual basis we no longer need a primal basis.
The programs BasisManipulationVerbose
and BasisManipulation
in the examples illustrate how to use these functions and make speed comparisons. 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 class 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_FPZZflex.h
module. These optional parameters are allowed to take non-default values only when Int==ZZ
and precision = DOUBLE
. 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 method.
long LatticeTester::LLLConstruction0 | ( | NTL::matrix< long > & | gen, |
const double | delta, | ||
long | r, | ||
long | c, | ||
NTL::vector< double > * | sqlen ) |
long LatticeTester::LLLConstruction0 | ( | NTL::matrix< NTL::ZZ > & | gen, |
const double | delta, | ||
long | r, | ||
long | c, | ||
NTL::vector< double > * | sqlen ) |
long LatticeTester::LLLConstruction0 | ( | NTL::matrix< NTL::ZZ > & | gen, |
const double | delta, | ||
long | r, | ||
long | c, | ||
NTL::vector< NTL::RR > * | sqlen ) |
long LatticeTester::LLLConstruction0 | ( | NTL::matrix< NTL::ZZ > & | gen, |
const double | delta, | ||
long | r, | ||
long | c, | ||
NTL::vector< quad_float > * | sqlen ) |
long LatticeTester::LLLConstruction0 | ( | NTL::matrix< NTL::ZZ > & | gen, |
const double | delta, | ||
long | r, | ||
long | c, | ||
NTL::vector< xdouble > * | sqlen ) |
void LatticeTester::lowerTriangularBasis | ( | IntMat & | gen, |
IntMat & | basis, | ||
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 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.
|
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 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.
|
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.
void LatticeTester::mDualBasis | ( | const IntMat & | basis, |
IntMat & | basisDual, | ||
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.
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.
void LatticeTester::mDualBasis | ( | const NTL::matrix< NTL::ZZ > & | basis, |
NTL::matrix< NTL::ZZ > & | basisDual, | ||
const NTL::ZZ & | m ) |
void LatticeTester::mDualUpperTriangular | ( | const IntMat & | A, |
IntMat & | B, | ||
const Int & | m, | ||
long | dim ) |
For B
to be m
-dual to A
, we have to have that \(AB^t = mI\).
Takes an upper triangular basis matrix basis
and computes the m-dual basis basisDual
.
Since 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}.\]
|
static |
Takes an upper triangular basis matrix basis
and computes the m-dual basis basisDual
.
The method assumes that each coefficient on the diagonal of basis
is nonzero and divides m
. That is, the basis matrix must be square and invertible. The algorithm is described in the Lattice Tester guide [7]. Since the basis is upper triangular, its m-dual will be lower triangular. When dim > 0
, it gives the number of rows and columns of the matrix basis
that is actually used. Otherwise (by default) the method uses basis.numCols()
.
Takes an upper triangular basis matrix basis
and computes the m-dual basis basisDual
.
Since 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}.\]
void LatticeTester::mDualUpperTriangular96 | ( | IntMat & | basis, |
IntMat & | basisDual, | ||
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].
|
static |
This function does essentially the same thing as mDualUpperTriangular
, but the algorithm is slightly different.
It uses the method described in [3].
void LatticeTester::mDualUpperTriangular96 | ( | NTL::matrix< NTL::ZZ > & | basis, |
NTL::matrix< NTL::ZZ > & | basisDual, | ||
const NTL::ZZ & | m, | ||
long | dim ) |
|
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 wierd type convertion and might not work well if different types are used.
|
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.
void LatticeTester::MyExit | ( | int64_t | status, |
std::string | msg ) |
Special exit function.
status
is the status code to return to the system, msg
is the message to print upon exit.
|
inline |
Returns \(\sqrt{x}\) for \(x\ge0\), and \(-1\) for \(x < 0\).
void LatticeTester::negativeCholesky | ( | ) |
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 }
.
std::ostream & LatticeTester::operator<< | ( | std::ostream & | out, |
const std::pair< T1, T2 > & | x ) |
Streaming operator for vectors.
Formats a pair as: (first,second)
.
std::ostream & LatticeTester::operator<< | ( | std::ostream & | out, |
const std::set< K, C, A > & | x ) |
Streaming operator for sets.
Formats a set as: { val1, ..., valN }
.
std::ostream & LatticeTester::operator<< | ( | std::ostream & | out, |
const std::vector< T, A > & | x ) |
Streaming operator for vectors.
Formats a vector as: [ val1, ..., valN ]
.
|
related |
void LatticeTester::printBase | ( | IntMat | bas_mat | ) |
void LatticeTester::printBase2 | ( | IntMat | bas_mat | ) |
|
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 many places to compute a floating point norm of vectors with integers values. Take care when using this function.
void LatticeTester::ProductDiagonal | ( | const NTL::matrix< 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
.
void LatticeTester::projectionConstruction | ( | const IntMat & | inBasis, |
IntMat & | projBasis, | ||
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.
|
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.
void LatticeTester::projectionConstructionLLL | ( | const IntMat & | inBasis, |
IntMat & | projBasis, | ||
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.
|
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.
void LatticeTester::projectionConstructionUpperTri | ( | const IntMat & | inBasis, |
IntMat & | projBasis, | ||
const Coordinates & | proj, | ||
const Int & | m, | ||
long | r ) |
|
static |
void LatticeTester::projectionConstructionUpperTri | ( | const IntMat & | inBasis, |
IntMat & | projBasis, | ||
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 method 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.
|
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 method 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.
void LatticeTester::projectMatrix | ( | const IntMat & | in, |
IntMat & | out, | ||
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.
|
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.
|
inline |
Computes a/b
, truncates the fractionnal 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\) |
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.
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.
std::int64_t LatticeTester::RandInt | ( | std::int64_t | i, |
std::int64_t | j ) |
double LatticeTester::RandU01 | ( | ) |
Returns a random number in \([0, 1)\).
The number will have 53 pseudo-random bits.
void LatticeTester::redBKZ | ( | IntMat & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
RealVec * | sqlen ) |
|
static |
This calls the NTL implementation of the floating point version of the BKZ reduction algorithm presented in [10], 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.
void LatticeTester::redBKZ | ( | NTL::matrix< int64_t > & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
NTL::vector< double > * | sqlen ) |
void LatticeTester::redBKZ | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
NTL::vector< double > * | sqlen ) |
void LatticeTester::redBKZ | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
NTL::vector< NTL::RR > * | sqlen ) |
void LatticeTester::redBKZ | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
NTL::vector< quad_float > * | sqlen ) |
void LatticeTester::redBKZ | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | blocksize, | ||
long | prune, | ||
long | dim, | ||
NTL::vector< xdouble > * | sqlen ) |
void LatticeTester::redLLL | ( | IntMat & | basis, |
double | delta = 0.99999, | ||
long | dim = 0, | ||
RealVec * | sqlen = 0 ) |
This file provides only static functions to reduce a lattice basis in some ways (pairwise, LLL, BKZ [4], [9], [10]).
The LLL and BKZ ones are wrappers to slightly modified versions of the NTL functions described at https://libntl.org/doc/LLL.cpp.html . These functions do not require the creation of an object like in ReducerBB
. They take an IntMat
object that contains a basis and return the reduced basis in the same object. In our versions, the basis can occupy only part of the IntMat
object so we can use the same object for several bases of various sizes, the basis entries can be of either ZZ or int64_t
type, the shortest basis vector is always placed in the first row, and the squared vector lengths are also returned in a vector. The norm to measure the vector lengths is always taken as the Euclidean one. The computations can be done either in double
or in RR
for the real numbers.
All these functions take as input (and return as output) an IntMat
object whose first dim
rows and columns (a square matrix) are assumed to form be a set of independent vectors that form a basis for the lattice. The same IntMat
object can then be used for several lattices of different sizes. This function uses the NTL implementation of the LLL reduction algorithm with factor delta
, presented in [10] (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 between 1/2 and 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.
|
static |
This file provides only static functions to reduce a lattice basis in some ways (pairwise, LLL, BKZ [4], [9], [10]).
The LLL and BKZ ones are wrappers to slightly modified versions of the NTL functions described at https://libntl.org/doc/LLL.cpp.html . These functions do not require the creation of an object like in ReducerBB
. They take an IntMat
object that contains a basis and return the reduced basis in the same object. In our versions, the basis can occupy only part of the IntMat
object so we can use the same object for several bases of various sizes, the basis entries can be of either ZZ or int64_t
type, the shortest basis vector is always placed in the first row, and the squared vector lengths are also returned in a vector. The norm to measure the vector lengths is always taken as the Euclidean one. The computations can be done either in double
or in RR
for the real numbers.
All these functions take as input (and return as output) an IntMat
object whose first dim
rows and columns (a square matrix) are assumed to form be a set of independent vectors that form a basis for the lattice. The same IntMat
object can then be used for several lattices of different sizes. This function uses the NTL implementation of the LLL reduction algorithm with factor delta
, presented in [10] (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 between 1/2 and 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.
void LatticeTester::redLLL | ( | NTL::matrix< int64_t > & | basis, |
double | delta, | ||
long | dim, | ||
NTL::vector< double > * | sqlen ) |
void LatticeTester::redLLL | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | dim, | ||
NTL::vector< double > * | sqlen ) |
void LatticeTester::redLLL | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | dim, | ||
NTL::vector< NTL::RR > * | sqlen ) |
void LatticeTester::redLLL | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | dim, | ||
NTL::vector< quad_float > * | sqlen ) |
void LatticeTester::redLLL | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta, | ||
long | dim, | ||
NTL::vector< xdouble > * | sqlen ) |
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).
|
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).
void LatticeTester::redLLLExact | ( | NTL::matrix< NTL::ZZ > & | basis, |
double | delta ) |
|
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.)
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.
|
inline |
Sets the first d
components of A
to 0.
Sets the first d
components of A
to the value x
.
|
inline |
Sets the first d
of A
to 0.
|
inline |
Sets the first d
components of A
to 0.
|
inline |
Returns the sign of x
.
The sign is 1 if x>0
, 0 if x==0
and -1 if x<0
.
|
static |
Helper function to skip all characters of a given class.
|
inline |
|
inline |
|
inline |
|
inline |
|
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.
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
.
std::string LatticeTester::toString | ( | Chrono & | timer | ) |
Returns the value of the duration from timer
.
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.
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*)
.
|
static |
|
static |
|
static |
|
static |
|
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.
|
static |
|
static |
|
static |
|
static |
|
static |
|
static |
|
static |
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
void LatticeTester::Triangularization | ( | Matr & | W, |
Matr & | V, | ||
int64_t | lin, | ||
int64_t | col, | ||
const Int & | m ) |
Takes a set of generating vectors in the matrix W
and iteratively transforms it into an upper triangular lattice basis into the matrix V
.
W
and V
have to have more rows than lin
and more columns than col
since this algorithm will only operate on the upper lin*col
matrix of W
. All the computations will be done modulo m
, which means that you must know the rescaling factor for the vector system to call this function. After the execution, W
will be a matrix containing irrelevant information and V
will contain an upper triangular basis.
For more details please look at [3]. This algorithm basically implements what is written at the end of the article, that is the matrix W
contains the set of vectors that is used and modified at each step to get a new vector from the basis.
void LatticeTester::upperTriangularBasis | ( | IntMat & | gen, |
IntMat & | basis, | ||
const Int & | m, | ||
long | dim1, | ||
long | dim2 ) |
Same as lowerTriangularBasis
, except that the returned basis is upper triangular.
|
static |
Same as lowerTriangularBasis
, except that the returned basis is upper triangular.
const double LatticeTester::MAX_LONG_DOUBLE = 9007199254740992.0 |
Maximum integer that can be represented exactly as a double
: \(2^{53}\).
const std::int64_t LatticeTester::TWO_EXP |
Table of powers of 2: TWO_EXP[
\(i\)]
\(= 2^i\), \(i=0, 1, …, 63\).