Magma provides an experimental implementation of the fastest general purpose factoring algorithm known: the number field sieve (NFS). In order to make use of the NFS the user needs to have some knowledge of the algorithm. The implementation requires a significant amount of memory and disk space. For factorizations of 80-digits or more, at least 64 megabytes of RAM and half a gigabyte of disk space is required.
The implementation can be used for both general number field and special number field factorizations: the only difference is in the polynomial selection. Presently, Magma does not provide a function to choose a good polynomial for a particular number to be factored. However, Magma does provide some functions that are useful for the implementation of the polynomial selection algorithms developed by Peter Montgomery and Brian Murphy in [Mur99]. These are described in the next section.
Magma's NFS does line-by-line sieving only (also known as "classical sieving"). Future versions will include lattice sieving. Also, Magma's NFS is restricted to one linear polynomial (the "rational side") and one polynomial of higher degree (the "algebraic side"). At the time of writing, this is not a major restriction since the best methods for selecting polynomials for numbers of at least 100-digits are of this form.
To factor a number n with Magma's NFS, the user must at the very least specify n, an algebraic polynomial F in homogeneous form, and two integers m1 and m2 such that F(m1, m2) = 0 mod n. The user can then call NFS(n, F, m1, m2) to perform a complete NFS factorization. However, Magma's default parameter choices are not very optimized. Almost certainly, the user can factor much faster by choosing his or her own parameters to optimize for their specific polynomial. Therefore, one should use this form of the implementation only to get a rough idea on parameter choices, and no conclusions should be drawn about the speed of the NFS when it is called in this way.
More precise control over the behaviour of NFS is achieved by calling the following six functions in the order: NFSRelations, NFSCycleFile, NFSCharacterColumns, NFSDependencies, NFSFactor, and NFSClear. Each of these functions requires the following input parameters: n, F, and a tuple containing the following data which is explained in more detail below: m1, m2, upper bound for algebraic primes in factor base, upper bound for rational primes in factor base, output data file name, an approximate value for the maximum value of a sieved out to, the first value of b sieved, the last value of b sieved, an upper bound for large primes on the algebraic side, an upper bound for large primes on the rational side, an error term for the sieve threshold on the algebraic side, and an error term for the sieve threshold on the rational side. An optional 13th element can be give in the tuple to take advantage of the cache memory size of the computer. This value should be 1 to indicate small cache size, 2 to indicate medium cache size, or 3 to indicate large. Most of the parameters from the tuple are only relevant for the sieving stage. The next few paragraphs explain them in more detail.
The line-by-line siever sieves values of F(a, b) on the algebraic side and corresponding values a ⋅m2 - b ⋅m1 on the rational side. This is done by fixing a value of b (starting with the first value of b specified in the input parameters), then sieving all values of a such that the -max_a <= a < max_a, where max_a is approximately equal to the sixth element of the tuple in the input parameters (some rounding off is done to make sure that the sieve interval length is divisible by a high power of 2). After this value of b is completed, b is incremented and the next value of b is processed. This continues until the maximum value of b is completed (specified in the input parameters), or until enough relations are obtained to complete the factorization. If the first value of b specified in the input parameters is 0, then free relations will be output to the data file and then sieving will start at b=1. If the last value of b specified in the input parameters is 0, then sieving will continue until enough relations are obtained to complete the factorization. Cycles are counted after every 256 iterations.
For the sieve, the rounded off natural logarithms of primes are used to mark the sieve interval. Moreover, the implementation does not sieve with prime powers. Hence, the code must allow for some error in scanning the sieve arrays for useful relations (i.e. the sieve threshold error input parameters). If, in addition, the user wants to take advantage of large prime relations (recommended), then larger error terms should be used. The implementation will keep relations having up to 2 large primes on each side, but will only find such relations if the user selects large enough sieve threshold error bounds. Be cautious, however, since using such relations greatly increases the disk space requirements. It takes a bit of experimentation to determine the best error bounds for speed or disk space optimization purposes.
During the sieve stage, relations that fully factor into primes within the factor bases ("full relations") will be output to the data file specified by the user in the tuple. Other useful relations ("partial relations") will be output to a separate file which has the name of the specified input file concatenated with "_partials". Whenever cycles [LD95] are counted, another file is created which has the name of the specified input file concatenated with "_cycles". Some other files are created and then deleted during the cycle counting process.
If any of the following input parameters are set to 0, default values will be chosen: the upper bound for algebraic primes in factor base, the upper bound for rational primes in factor base, the approximate value for the maximum value of a sieved out to, the upper bound for large primes on the algebraic side, or the upper bound for large primes on the rational side. Using the default factor base bounds could give very poor performance, so it is recommended that these values are used only as a starting point in determining good parameters for a particular number.
For users wanting to distribute the computations among several processors, this is only possible for the sieving stage and the final factorization stage. To do this for the sieving stage, each processor should get a unique range of b-values to sieve and unique data file names. During the sieving, the user must manually check when the combined data has enough relations to factor the number. To do this, the data files must first be merged using NFSMerge and then the cycles can be counted with NFSCycleCount. If the combined number of full relations plus the number of cycles exceeds the size of both factor bases combined, then sieving can discontinue, and the user can proceed to the other stages of the factorization attempt using the merged data file name. Parallelizing the final factorization stage can be done by simply specifying in the NFSFactor command which matrix dependency each processor should use.
At the time of writing, the record NFS factorizations were lead by the CWI group and by people using CWI's or Arjen Lenstra's code. For those wanting to use Magma to assist in CWI factorizations, the command NFSCWIFormat can be used to convert data files from the Magma format to the CWI format. The resulting data file will be the input file name concatenated with "_CWI_format", and will contain both the full and partial relations.
The functions BaseMPolynomial, MurphyAlphaApproximation, OptimalSkewness, BestTranslation, PolynomialSieve, and DickmanRho, will be useful for those wanting to implement their own polynomial selection routines within the Magma interpreter language. These functions are described in the next section.
As a guideline for the selection of parameters, we include here a few examples of parameters that we have determined to be good for the Magma NFS implementation.
The best choice of factor base sizes depends on many variables, including the average log size and the alpha parameter for the polynomial F. Our polynomials above are quite good: if the user does not know much about the quality of polynomials, then he or she should use much larger factor bases.
> n := 5235869680233366295366904510725458053043111241035678897933802235060927; > P<X,Y> := PolynomialRing(Integers(),2); > F := 2379600*X^4 - 12052850016*X^3*Y - 13804671642407*X^2*Y^2 + > 11449640164912254*X*Y^3 + 7965530070546332840*Y^4 ; > m1 := 6848906180202117; > tuple := <m1,1,8*10^5,6*10^5,"/tmp/nfs_data",4194280,0,0,0,0,16,14>;
> n := 1871831866357686493451122722951040222063279350383738650253906933489072\ > 2483083589; > P<X,Y> := PolynomialRing(Integers(),2); > F := 15901200*X^4 + 9933631795*X^3*Y - 112425819157429*X^2*Y^2 - > 231659214929438137*X*Y^3 - 73799500175565303965*Y^4; > m1 := 1041619817688573426; > tuple := <m1,1,800000,600000,"/tmp/nfs_data",10485760,0,0,0,0,20,16>;
> n := 12118618732463427472219179104631767765107839384219612469780841876821498\ > 2402918637227743; > P<X,Y> := PolynomialRing(Integers(),2); > F := 190512000*X^4 - 450872401242*X^3*Y + > 1869594915648551*X^2*Y^2 + 2568544235742498*X*Y^3 - > 9322965583419801010104*Y^4; > m1 := 28241170741195273211; > tuple := <m1,1,1600000,1000000,"/tmp/nfs_data",2^24,0,0,0,0,24,18>;
Performs the factorization of n using the number field sieve with algebraic polynomial F, where m1 and m2 are such that F(m1, m2) = 0 mod n. Returns a nontrivial factorization of n if one is found, or 0 if no factor is found. Note that the default algebraic factor base size of NFS is chosen to be rather large to decrease the likelihood of running out of useful relations. Unless the user chooses a really bad polynomial, this should not happen.
Generates the relations to factor n by the number field sieve factoring algorithm using algebraic homogeneous polynomial F. tuple should be of the form < m1, m2, y_alg, y_rat, data_file, max_a, first_b, last_b, lp_bound_alg, lp_bound_rat, alg_error_term, rat_error_term > where m1 and m2 are such that F(m1, m2) = 0 mod n, y_alg is upperbound for algebraic primes in factor base, y_rat is upperbound for rational primes in factor base, data_file is the name of a file that the full relations will be put into, max_a is the approximate maximum value of a sieved for each value of b, first_b is the first value of b sieved, last_b is the last value of b sieved, lp_bound_alg is an upperbound for large primes on the algebraic side, lp_bound_rat is an upperbound for large primes on the rational side, and alg_error_term and rat_error_term are error terms used in scanning the sieve array for useful relations on the algebraic side and rational side respectively. An optional 13th element can be given for tuple to take advantage of the cache memory size of the computer: a value of 1 indicates small cache size, 2 for medium cache size, and 3 for large cache size. Returns the number of full relations plus the number of cycles found.
Merges the NFS data files listed in T into new data files, while removing duplicate and corrupted lines of data. The full relations go into a file with name fn, and the partial relations into a file named fn with "_partials" appended to it. This function is mainly intended for factoring with multiple processors or to clean up data files that may have errors in them due to a process being interrupted.
Converts data files to CWI format, storing primes only greater than or equal to pb. The resulting data file name will be the same as the input file name (specified in T) with "_CWI_format" appended to it.
Returns the number of cycles in the partial data file corresponding to the file name specified in tuple. This function is mainly intended for factoring with multiple processors.
Creates a file with all the cycle information that the Magma NFS needs for the matrix and final factorization stages. The input parameters are the same as in NFSRelations.
Creates a file with the quadratic character column data for the full relations and the cycles. If c is specified, then 32c quadratic character columns are created. Otherwise, 32 character columns are created. The other input parameters are the same as in NFSRelations.
Uses block Lanczos to find dependencies in the NFS matrix. Dependencies will be output to a data file. The input parameters are the same as in NFSRelations.
Attempt to factor n using the output from NFSDependencies. If k is specified, then use the k'th dependency output from the matrix, where the first dependency is number 0. Otherwise, attempt to use all dependencies until a nontrivial factorization is found. The other input parameters are the same as in NFSRelations. Returns a nontrivial factor if successful, or else 0 if unsuccessful.
Deletes the data files that were created during various stages of the NFS implementation.
Returns a homogeneous degree d polynomial in two variables which is obtained by writing n in base m and adjusting coefficients to be between -m/2 and +m/2. The return polynomial with m1 = m and m2 = 1 can be used to obtain the factorization of n with the number field sieve. Must have d >= 2 and n >= m^d.
Approximates the alpha (as described in Brian Murphy's PhD Thesis) value for the polynomial F using primes up to b. F may be in univariate or homogeneous form. Since random sampling is used for primes dividing the discriminant, multiple calls to this function will give slightly different answers.
Return the optimal skewness value for the input polynomial, and the corresponding average log size. These terms are explained in Brian Murphy's PhD thesis. F may be in univariate or homogeneous form.
> n := RandomPrime(90)*RandomPrime(90);
> n;
3596354707256253204076739374167770148715218949803889
> d := 4;
> approx_m := Iroot( Iroot( n, d+1 ) * Iroot( n, d ) , 2 );
> leading_coeff := n div approx_m^d;
> leading_coeff;
143082
> m := Iroot( n div leading_coeff, d );
> P<X,Y> := PolynomialRing( Integers(), 2 );
> F<X,Y> := BaseMPolynomial(n,m,d);
> F;
143082*X^4 + 463535*X^3*Y - 173869838910*X^2*Y^2 + 167201617413*X*Y^3 +
159859288415*Y^4
> skew, als := OptimalSkewness( F );
> alpha := MurphyAlphaApproximation( F, 2000 );
> rating := als + alpha;
> rating;
23.143714548914575193314917
>
> best_rating := rating;
> best_m := m;
> for i in [1..100] do
> leading_coeff := leading_coeff + 1;
> m := Iroot( n div leading_coeff, d );
> F<X,Y> := BaseMPolynomial(n,m,d);
> skew, als := OptimalSkewness( F );
> alpha := MurphyAlphaApproximation( F, 2000 );
> rating := als + alpha;
> if rating lt best_rating then
> best_rating := rating;
> best_m := m;
> end if;
> end for;
> best_rating;
20.899568473033257031950385
> best_m;
398116527578
> F<X,Y> := BaseMPolynomial(n,best_m,d);
> F;
143160*X^4 + 199085*X^3*Y - 9094377652*X^2*Y^2 - 93898749030*X*Y^3 -
169859083883*Y^4
> OptimalSkewness( F );
165.514255523681640625 20.969934467920612180646408
> MurphyAlphaApproximation( F, 2000 );
-0.0542716157630141449500150842
> time NFS( n, F, best_m, 1 );
...
Given a tuple that contains a univariate polynomial f, an integer m, a real number corresponding to the average log size of f for some optimal skewness value s, and optionally the real number s, this routine finds a translation of f such that the average log size is a local minimum. Returns a tuple containing the translated polynomial g, an integer m' such that g(m') = f(m), the average log size of g over a rectangle having optimal skewness s', and s'. This function also accepts the homogeneous forms of polynomials, and when called in this way, the return polynomial will also be in homogeneous form.
Given a tuple containing a univariate polynomial f, an integer m, an integer J0, an integer J1, and a real number gamma, this function searches for polynomials of the form g = f + (j_1 x - j_0)(x - m) with | j_0 | < J0, | j_1 | < J1, and having alpha value no larger than gamma. The search is done via a sieve as described in Brian Murphy's PhD thesis. Returns a sequence of tuples, where each tuple contains a rating of the form the approximate alpha plus the average log size for some optimal skewness s of some polynomial g, the approximate alpha, the average log size, s, g, m, j_0, and j_1. Smaller ratings tend to indicate better polynomials.
Computes rho(u) where rho is Dickman's rho function.[Next][Prev] [Right] [Left] [Up] [Index] [Root]