Monday 31 December 2007

Arjen Bot

Also there is Arjen Bot's page of factorizations of 2^n+1
http://www.euronet.nl/users/bota/medium-p-odd.txt

Sunday 30 December 2007

Saturday 29 December 2007

NFSNET

Here is a link to the NFSNET project, a distributed attack on Cunningham Project targets using the Number Field Sieve:
http://www.nfsnet.org/

Friday 28 December 2007

P49_101_87

I can't resist mentioning a fun result I had recently for the XYYXF project:
http://xyyxf.at.tut.by/records.html#ecm
[There I am currently in tenth place! :)]
I can certainly recommend the XYYXF ECMNet server as a slick and easy way to find interesting factors for that project.

Thursday 27 December 2007

WE on Carmichaels

When running WE on Carmichaels, clearly the algorithm will never terminate.
http://en.wikipedia.org/wiki/Carmichael_numbers
However, there is a simple workaround, by simply multiplying the input number to be factored (if Carmichael) by a suitable larg(ish) prime, eg a Mersenne. This will normally eject the factors, additionally, as required.
Here are some examples:
1)
tiggatoo:~/math/we james$ java we2tr2
f168003672409*(2^127-1)

[28584343649372241809274301558307881708368828786343]
base#72
elapsed=0s factor=6073 A=13044086312830013543739988769
base#70
elapsed=0s factor=3037 A=860763419187377590694240501658
base#230
elapsed=1s factor=9109 A=771019052198021707484273222802
170141183460469231731687303715884105727
duration: 1 seconds


2)
tiggatoo:~/math/we james$ java we2tr2
f173032371289*(2^89-1)

[107101850255573582977951074701018631079]
base#43
elapsed=0s factor=3067 A=1236330343540696294158885383112
base#149
elapsed=0s factor=6133 A=559525871748815289923955370298
base#758
elapsed=1s factor=9199 A=194689393018084379227390724801
618970019642690137449562111
duration: 1 seconds


3)
tiggatoo:~/math/we james$ java we2tr2
f973694665856161*(2^127-1)

[165665562777913375106491436985024163141370898298334047]
base#3
elapsed=0s factor=6841 A=1245204978926759694403655173839
base#0
elapsed=0s factor=2281 A=1024632434440788626941999868632
base#7
elapsed=0s factor=4561 A=457244767350940464517786915322
base#12
elapsed=0s factor=13681 A=85779304638188382542657979233
170141183460469231731687303715884105727
duration: 0 seconds


And here, for the record/or in case anyone wants to confirm results above, is the Java source code to we2tr2 (Copyright JGW ~ 2000):

import java.math.BigInteger;
import java.util.Random;
import java.util.Date;

public class we2tr2 {
static BigInteger zero = new BigInteger("0");
static BigInteger one = new BigInteger("1");
static BigInteger two = new BigInteger("2");
static BigInteger hundred = new BigInteger("100");
static BigInteger thousand = new BigInteger("1000");

static BigInteger n = new BigInteger("0");
static BigInteger p = new BigInteger("0");
static BigInteger known = new BigInteger("0");
static BigInteger numtrials = new BigInteger("0");
static String s;

static int olds1 = 0;
static int s1 = 0;

static Date d;
static long starttime;
static long finishtime;
static long duration;




public static void main (String args[])
throws java.io.IOException {

we2tr2 we2tr2inst = new we2tr2();

char c;
String sInput;

StringBuffer sbInput = new StringBuffer("");

while ((c = (char)System.in.read()) != '\n' && c != '\r')
sbInput.append(c);
System.in.read();
sInput = sbInput.toString().trim();


if (sInput.charAt(0) == 'f' || sInput.charAt(0) == 'F') {
s = sInput.substring(1).trim();

s1 = 0;
olds1 = 0;
p = we2tr2inst.eval(s);

System.out.println('[' + p.toString() + ']');
}
else {
p = new BigInteger(sInput);
}


n = p;

d = new Date();
starttime = d.getTime();

we2tr2inst.factorize(n);

d = new Date();
finishtime = d.getTime();

duration = (finishtime-starttime)/1000;

System.out.println("duration: " + duration + " seconds");

System.out.println();
}

public boolean factorize(BigInteger n) {
boolean prime = false;
BigInteger numtested = new BigInteger("0");
BigInteger T = new BigInteger("1");
BigInteger b = new BigInteger("1");
BigInteger A = new BigInteger("2");
BigInteger wanless = new BigInteger("2");

if (n.isProbablePrime(1000)) {
prime = true;
System.out.println(n);
return prime;
}

// workaround - apparent java bug in modPow - JGW
if (n.compareTo(two) < 0)
return false;

if (n.remainder(two).compareTo(zero) == 0) {
System.out.println(two.toString());
// added 2006-06-09
return(factorize(n.divide(two)));
}
// end workaround

while (wanless.compareTo(n) < 0)
wanless = wanless.multiply(two);

Random r = new Random();

numtested = zero;

while (T.compareTo(one) == 0 || T.compareTo(n) == 0) {
// changed JW 2005-3-23
A = new BigInteger(hundred.intValue(), r);

// added JGW 2006-06-09
System.out.print("base#" + numtested + '\r');



// changed DT 2005-2-20
b = A.modPow(wanless, n);
T = n.gcd(b.modPow(n, n).subtract(b));

numtested = numtested.add(one);
}





if (T.compareTo(one) > 0 && T.compareTo(n) < 0) {

d = new Date();
finishtime = d.getTime();

duration = (finishtime-starttime)/1000;


System.out.println();
System.out.println("elapsed=" + duration + "s" + '\t' + "factor=" + T.toString() + '\t' + "A=" + A.toString() + '\t');
factorize(n.divide(T));
}

return prime;

}


public BigInteger evalRand(char op, BigInteger oldn) {
BigInteger n = new BigInteger("1");


switch (op) {
case 'r':
case 'R':
Random r = new Random();

n = new BigInteger(oldn.intValue(), r);
break;

default:
n = oldn;
break;
}
return n;
}

public BigInteger evalFact(BigInteger oldn, char op) {
BigInteger n = new BigInteger("1");
BigInteger i = new BigInteger("1");
BigInteger j = new BigInteger("1");
boolean prime = true;


switch (op) {
case '!':
for (i = one; i.compareTo(oldn) <= 0; i = i.add(one))
n = n.multiply(i);
break;

case '#':
for (i = one; i.compareTo(oldn) <= 0; i = i.add(one)) {
prime = true;
for (j = two; (prime == true) && (j.multiply(j).compareTo(i) <= 0); j = j.add(one))
if (i.remainder(j).compareTo(zero) == 0)
prime = false;
if (prime == true)
n = n.multiply(i);
}
break;

default:
n = oldn;
break;
}
return n;
}

public BigInteger evalPower(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '^':
n = oldn.pow(n1.intValue());
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger evalProduct(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '*':
n = oldn.multiply(n1);
break;
case '/':
n = oldn.divide(n1);
break;
case '%':
n = oldn.remainder(n1);
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger evalSum(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '+':
n = oldn.add(n1);
break;
case '-':
n = oldn.subtract(n1);
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger eval(String s) {
BigInteger oldn0 = new BigInteger("0");
BigInteger oldn1 = new BigInteger("0");
BigInteger oldn2 = new BigInteger("0");
BigInteger n = new BigInteger("0");

char oldop0 = 0;
char oldop1 = 0;
char oldop2 = 0;
char op = 0;

while (s1 < s.length()) {
switch (s.charAt(s1)) {
case '(':
case '[':
case '{':
s1++;
n = eval(s);
break;

case '0':
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
n = readNum(s);
break;

default:
break;
}
if (s1 < s.length()) {
switch (s.charAt(s1)) {

case ')':
case ']':
case '}':
case '!':
case '#':
case 'r':
case 'R':
case '^':
case '*':
case '/':
case '%':
case '+':
case '-':
op = s.charAt(s1);
s1++;
break;

default:
break;
}
}
else
op = 0;


switch (op) {
case 0:
case ')':
case ']':
case '}':
n = evalPower(oldn2, n, oldop2);
n = evalProduct(oldn1, n, oldop1);
n = evalSum(oldn0, n, oldop0);
return n;

case '!':
case '#':
n = evalFact(n, op);
break;

case 'r':
case 'R':
n = readNum(s);
n = evalRand(op, n);
break;

case '^':
n = evalPower(oldn2, n, oldop2);
oldn2 = n;
oldop2 = op;
break;

case '*':
case '/':
case '%':
n = evalPower(oldn2, n, oldop2);
oldop2 = 0;
n = evalProduct(oldn1, n, oldop1);
oldn1 = n;
oldop1 = op;
break;

case '+':
case '-':
n = evalPower(oldn2, n, oldop2);
oldop2 = 0;
n = evalProduct(oldn1, n, oldop1);
oldop1 = 0;
n = evalSum(oldn0, n, oldop0);
oldn0 = n;
oldop0 = op;
break;
default:
break;
}
}
return n;
}

public BigInteger readNum(String s) {
BigInteger n = new BigInteger("0");

olds1 = s1;
while (s1 < s.length() && Character.isDigit(s.charAt(s1)))
s1++;
n = new BigInteger(s.substring(olds1, s1));

return n;
}


}

Wednesday 26 December 2007

Charles Babbage

Here is a sketch (from Wikipedia) of Charles Babbage (b 26 December 1791).
He was the first to envisage a machine that could compute, though at the time this vision was of a purely mechanical device.
From http://en.wikipedia.org/wiki/Babbage
"He began in 1822 with what he called the difference engine, made to compute values of polynomial functions"
"The first difference engine was composed of around 25,000 parts, weighed fifteen tons (13,600 kg), and stood 8 ft (2.4 m) high. Although he received ample funding for the project, it was never completed"
"In 1991 a perfectly functioning difference engine was constructed from Babbage's original plans. Built to tolerances achievable in the 19th century, the success of the finished engine indicated that Babbage's machine would have worked"

Tuesday 25 December 2007

P2007

To celebrate Christmas 2007, for the past few days I have been factorizing P2007.
Here are the results:
ibu:~/math james$ time java superfac9
f2^2007+1

[14696072899510457910180264975074329395485666586735298566113827031369808145822340017365241424851280254956108347379039523500123122699047108242251921358933160773008638610599971840088163730974725743542902654728126239332046779346737710585256579333179693308275839559444787047544912589519783891140629020412202583212053620350010688717104574055412999539319651392054912347738448106306817040926244005345442289064602444671410741520258787821875717396461207456197233847539467765831034596299478021012490490523728714592688694474716929987628644661687302977141155300336976022455747686505323874664699578081559660947075760129]
wanless...
brutep: 3
wanless...
brutep: 3
wanless...
brutep: 3
brute...
brutep: 19
wanless...
ecm...
aprtcle: 247531
ecm...
aprtcle: 219256122131
wanless...
aprtcle: 20493495920905043950407650450918171260318303154708405513
ecm...
aprtcle: 4340301546362831119363
aprtcle:
56377694445208154141927654613855613062927113955212040908548454699046039020893338370875013074480485757794923
ecm...
aprtcle: 73215361
ecm... ^C
real 3459m6.204s
user 3424m7.740s
sys 5m57.800s

Monday 24 December 2007

Andrew Odlyzko

Andrew Odlyzko (picture reproduced with kind permission)
http://www.dtc.umn.edu/~odlyzko/
is probably mainly known for his contribution to the analysis of Riemann zeros, and the study of discrete logarithms. However, he also helped develop the Lanczos linear algebra stage of the QS (and other sieves), in his paper "Solving large sparse linear systems over finite fields" (1991),
http://citeseer.ist.psu.edu/140341.html
and has also authored a description of the state of factorization, looking ahead, "The future of integer factorization" (1995)
http://www.dtc.umn.edu/~odlyzko/doc/crypto.html

Saturday 22 December 2007

Tom Flowers

Today marks the birthday of Dr Tommy Flowers, the builder (somewhat w/ Alan Turing) of the first modern ie digital computer. It was constructed during WWII here in UK, to decrypt German operations messages, and christened 'Colossus'. It contained 2400 thermionic valves (picture from Wikipedia) each several centimetres long, so its name was suitably apt! These valves were basically operating as switches (even though they were primarily designed as amplifiers), and serve the same role that transistors on a chip do in modern microprocessors, though these latter contain many millions on a tiny wafer.
The Colossus represented a real advance in computing though, as it was the first time that a machine had been built using electronic rather than mechanical switches for greater speed and flexibility. It was this potential that Flowers recognised.

From http://www.ivorcatt.com/47c.htm
"At the time I had no thought or knowledge of computers in the modern sense and had never heard the term used except to describe somebody who did calculations on a desk machine."

"Colossus was useful in more than one way, and there were even demonstrations applying it to number theory. But these demonstrations were more notable for their ingenuity than for their effectiveness."

Friday 21 December 2007

Thursday 20 December 2007

Primorials(+/-1) Factorization

[The first in a trio of links, for this blog, of some types of numbers currently being factorized]
Factorization of Primorials(+/-1):
http://primorial.unit82.com/

Wednesday 19 December 2007

Factor Announcements

Here is a nice page detailing historical announcements of factorization achievements:
http://www.crypto-world.com/FactorAnnouncements.html

Monday 17 December 2007

Msieve v1.32

New version of msieve available:
Announcement at:
http://mersenneforum.org/showpost.php?p=120869&postcount=344
Download from:
http://www.boo.net/~jasonp/qs.html
It would appear that the promised merge w/ GGNFS is underway...

Sunday 16 December 2007

Number of University Departments Factorizing

How many University departments around the world are actively engaged in factorization? These spring to mind...
Well, there's Bruce Dodson's university, Lehigh in Pennsylvania USA for one - I understand Bruce manages several dozen PC's running ECM curves on the Cunningham Project targets.
Paul Zimmermann, at INRIA, in Nancy, France, is in charge of developing the GMP-ECM software itself, that Bruce Dodson (amongst others :) is running. I believe Alex Kruppa has links there, as well as many PC's in a French grid, the GRID5000. Amazingly, I'm guessing, the 5000 refers to an approximate, or at least, goal, of number of cores in the grid:
https://www.grid5000.fr/mediawiki/index.php/Special:G5KHardware
Greg Childers, at Cal State, Fullerton, California USA, also has several dozen machines available for factoring, which he uses on a variety of projects.
At Cambridge University, UK, Paul Leyland has, and still is, contributing to a variety of projects, including Cunningham Project (via NFSNET - more on that project at a later date, hopefully) and more recently, Homogeneous Cunninghams.
Then of course there is the (legendary?) CWI, or Centrum voor Wiskunde en Informatica in Amsterdam in the Netherlands...
There must be many more (maybe these are just the most vociferous!:), feel free to add to the list, you other folks out there [or let me know of specific omissions] - but that's 5 for starters, [I seem to have picked mainly computational rather than theoretical research] almost without thinking!

Saturday 15 December 2007

Tesla

...or you can buy a dedicated device (one example pictured (c) NVIDIA), called a 'Tesla', with the NVIDIA CUDA GPUs in situ:
http://en.wikipedia.org/wiki/NVIDIA_Tesla

Friday 14 December 2007

Msieve V1.31

New version of msieve (by Jason Papadopoulos) now available. It seems this is a significant upgrade, with lots of new stuff.
Download from:
http://www.boo.net/~jasonp/qs.html
Announcement here:
http://mersenneforum.org/showpost.php?p=120642&postcount=339

Thursday 13 December 2007

CUDA

Here are some links to an interesting new technology, which might be very useful for any "embarassingly parallel" (see earlier blog-post for definition of that) factorization (or other) algorithms. It's basically a method of using all the processing power inherent in modern GPUs for CPU-type calculations, and something I for one will be keeping an eye on...
http://en.wikipedia.org/wiki/CUDA
http://developer.nvidia.com/object/cuda.html#downloads
http://courses.ece.uiuc.edu/ece498/al1/
http://courses.ece.uiuc.edu/ece498/al1/Syllabus.html
It's even coming to Mac! :)
http://forums.nvidia.com/index.php?showtopic=47884

Wednesday 12 December 2007

GMP-M+2

I've recently noticed an interesting feature of GMP's implementation of mod_pow (ie modular exponentiation). I'm using this as part of random-base WEP on large M+2's. Anyway I always imagined, and hoped it would use "Russian Peasant" (see earlier post on this blog for description of that). If that is the case,
http://gmplib.org/manual/Powering-Algorithms.html#Powering-Algorithms
there should theoretically be a speed advantage to testing numbers that are (very) close to a power of 2 (as M+2's obviously are). And indeed I am observing just that - a 40% (!) speed increase testing the _raw_ M+2, rather than dividing out any factors (even though dividing out factors obviously leaves a slightly smaller number).

Tuesday 11 December 2007

C138_101_75

I'm a total xyyxf junkie! :)
Recently started another GNFS factorization for them, this time of 138-digit input number.
ETA several months... [on 1 CPU](assuming all goes well). Currently sieved 9502/839184 rels. required.
Hopefully the process will further increase my understanding of GNFS, which is kindof the point!

Monday 10 December 2007

WIFC

This entry also is probably somewhat overdue:
Check out Hisanori Mishima's comprehensive page on factorization at:
http://www.asahi-net.or.jp/%7EKC2H-MSM/mathland/matha1/
Note that there are also download links to factoring software (by Satoshi Tomabechi) on this page (especially MPQS-type factoring)

Sunday 9 December 2007

Cunningham Project

This is probably somewhat overdue, but here are a couple of links to the Cunningham Project, "one of the oldest continuously ongoing activities in computational number theory". [apparently it's been around since 1925]
http://en.wikipedia.org/wiki/Cunningham_project
http://homes.cerias.purdue.edu/~ssw/cun/index.html

Saturday 8 December 2007

So just how big are these numbers?

So just how big are these numbers?
Well mathematics is unique among the sciences for generating really SILLY (or impressive, depending on your point-of-view :) sized numbers. Even (standard) physics not only doesn't come close, but can't even stand comparison. Hence the other-worldly nature of math (and mathematicians! :)). Now possibly combinatorial aspects of holistic physics theories can generate big numbers, but in this case the theory is as much mathematical as physical, anyway.
The way math generates these numbers is by a trick called 'exponentiation' (usually). This allows one to write (and sometimes even calculate with and manipulate) a very long number in a very succinct form. Exponentiation is basically a shorthand notation for multiplying a number by itself many times. Thus 10^100 (also known as a googol) is the number you get by multiplying 10 by itself 100 times. If written out in full it is '1' followed by 100 zeroes. And this isn't even that big by math (or even factorization) standards. Some of the larger M+2 numbers I've been testing would have millions of digits if written out in full. You can see why I called these numbers 'silly'-sized! And if you want to go really crazy in math, you can even iterate/stack the exponents... (there are notations existent for this)
I expect many of you will be aware of the story of the Chinese inventor of chess? Apparently the emperor was so impressed with chess that he promised the inventor a reward based on the 64 squares of the chess-board, that the inventor (obviously a mathematician too) had tricked him into granting: namely this:
1 grain of rice for the first square, 2 for the second, 4 for the third, 8 for the fourth, and so on, doubling each time. Total number of grains of rice 2^64, or about 10^20 - far more than the whole production of China, and the emperor was never able to make good his promise...! This doubling process is another example of exponentiation, as the emperor clearly learned the hard way.
For comparison:
1) Number of seconds elapsed since Genesis ~ 2*10^11 (just a measly 12-digit number)
2) Number of seconds elapsed since the Jurassic ~ 2*10^15 (just a measly 16-digit number)
3) Number of seconds elapsed since the start of the visible universe ~ 3*10^20 (just a measly 21-digit number)
Incidentally, the complete factorization of any of these (exact) numbers on a modern PC would only take a (very) split-second.
I'll repeat...
Some of the numbers under test at the Mersenneplustwo project, for example, have _millions_ of digits. See why I used the adjective 'SILLY' yet - there really is almost no other word for these-size numbers?
Finally, and this is perhaps the most remarkable fact of all - mathematics also deals with the infinite. And compared to infinity (ie the _whole_ Universe) every single one of these numbers is actually infinitesimally SMALL! Now how do you get your head around that??? I know I for one, have problems with that...
Ahh, the mysteries of the finite and the infinite - or as Shakespeare famously put it:
"To be or not to be, that is the question"

Friday 7 December 2007

M+2 record

Mersenneplustwo project.
Here is a record of the number, and size, of new factors found, as well as an indication of the total effort, for each year since 2005:
2005 - 10 (largest 38-digits) - 4GHz-yrs
2006 - 3 (largest 19-digits) - 12GHz-yrs
2007[so far] - 2 (largest 23-digits) - 100GHz-yrs

Thursday 6 December 2007

M+2 progress

A bit of news from the Mersenneplustwo project:
As hoped for, getting Lenny the new MacBook going on those mid-range M+2's with mprime (ie by ECM) recently yielded a 21-digit factor from (M+2)1257787:
http://bearnol.is-a-geek.com/Mersenneplustwo/Mersenneplustwo.html
Thanks once again to George Woltman, for his cool software!
http://www.mersenne.org/freesoft.htm
This is the second new factor found this year (both by mprime)...

Wednesday 5 December 2007

John Tukey


...to conclude this series on FFTs:
This is a (PD) picture of John Tukey.

Tuesday 4 December 2007

FFTs-part5

The Cooley-Tukey Algorithm is an FFT algorithm especially suited, and adapted for large integer multiplication. It was first published in 1965.
GIMPS' mprime uses "radix-4" Cooley-Tukey to achieve a computation time of the order of NlogN (for N*N multiplication).
http://en.wikipedia.org/wiki/Cooley-Tukey_FFT_algorithm

Monday 3 December 2007

FFTs-part4

Parseval's Theorem
A special case of Plancherel's Theorem, when function x = function y. Then the scaling constant = 1/(2*pi) [FT] or 1/N [DFT]
http://en.wikipedia.org/wiki/Parseval%27s_theorem
"The interpretation of this form of the theorem is that the total energy contained in a waveform x(t) summed across all of time t is equal to the total energy of the waveform's Fourier Transform X(f) summed across all of its frequency components f."

Sunday 2 December 2007

FFTs-part3

Plancherel's Theorem
A (unitary) relationship between the DFTs of any two functions x and y under multiplication (subject to a constant scaling factor), which property can be of additional use in simplifying the problem under computation.
http://en.wikipedia.org/wiki/Plancherel_theorem

Saturday 1 December 2007

FFTs-part2

The DFT (discrete Fourier transform) is basically a digital sampling (ie approximation) of the analog Fourier transform, and as such it remains finite throughout the domain.
http://en.wikipedia.org/wiki/Discrete_Fourier_transform
A fast Fourier transform (FFT) is an efficient algorithm to compute the discrete Fourier transform and its inverse, for eg, a particular problem.
http://en.wikipedia.org/wiki/Fast_Fourier_transform

Friday 30 November 2007

FFTs-part1

The Fourier Transform - what is it?
Well, suppose we have some sum or computation we wish to evaluate, for example to calculate a*b for some a and b, e.g. large integers (of several thousands of digits or more).
Sometimes it can be easier to map each of the inputs to a new domain via some transform, perform the multiplication (approximating suitably if necessary) and then apply the reverse transform to get the original answer. Such a transform is the "Fourier" transform (and its inverse), and because of the way it maps objects/events to wave functions (expressed in terms of sin's and cos's) (and back) its respective domains of operation are termed 'time' and 'frequency'.
http://en.wikipedia.org/wiki/Fourier_transform
http://mathworld.wolfram.com/FourierTransform.html
http://mathworld.wolfram.com/FourierSeries.html
Obviously, in order for the process to work, it depends on the reciprocity of the transforming function. This is achieved by the orthogonality of the underlying individual waveforms (Sturm-Liouville).
http://en.wikipedia.org/wiki/Sturm-Liouville_theory

Thursday 29 November 2007

Plancherel


This is a (PD) picture of Michel Plancherel.
...which should lead me on to a discussion of Fourier Transforms.

Wednesday 28 November 2007

C127_113_36

Well, I finally completed my factorization of C127_113_36 for the XYYXF project - and by the recommended method - ie sieving by GGNFS, followed by post-processing with msieve. Thanks to Greg Childers, Bob Backstrom and Hallstein Hansen, who helped me with this transfer.
Basically msieve needs three files:
1) worktodo.ini - just containing the input number
2) msieve.fb - with the details (roughly) transferred from the .poly file of GGNFS
3) msieve.dat - with all the actual relations from GGNFS, translated to msieve format by procrels w/ the following crucial command:
"procrels -fb C127_113_36.fb -prel rels.bin -dump"

After this it's just a case of running msieve w/
"msieve -nc -v"
In the end the sieving [Lintel] took me about 3 months(!), because of memory limitations, and the GGNFS 'bounce', while the postprocessing in msieve [PPC-Tiger] finished to schedule in less than 2 days.

Tuesday 27 November 2007

Sieving Records

In fact, here is a nice page detailing historical records of factorizations using sieving methods:
http://www.crypto-world.com/FactorRecords.html

Monday 26 November 2007

Quadratic Sieve

Also, here is some info on the Quadratic Sieve (QS or MPQS, or, with slight variation, SIQS) method, invented in 1981 by Carl Pomerance, as an improvement to Dixon's Method:
http://en.wikipedia.org/wiki/Quadratic_sieve
http://mathworld.wolfram.com/QuadraticSieve.html
From the former link:
"On April 2, 1994, the factorization of RSA-129 was completed using QS. It was a 129-digit number, the product of two large primes, one of 64 digits and the other of 65. The factor base for this factorization contained 524339 primes. The data collection phase took 5000 MIPS-years, done in distributed fashion over the Internet. The data collected totaled 2GB. The data processing phase took 45 hours on Bellcore's MasPar (massively parallel) supercomputer. This was the largest published factorization by a general-purpose algorithm, until NFS was used to factor RSA-130, completed April 10, 1996."

Dixon's Method

And here is some info about Dixon's method, which is related to CFRAC, and the precursor to most other modern sieving methods.
http://en.wikipedia.org/wiki/Dixon%27s_factorization_method
http://mathworld.wolfram.com/DixonsFactorizationMethod.html
This method was first published in 1981.

Sunday 25 November 2007

CFRAC

Here are some links describing the CFRAC, or "continued fraction" factorization method:
http://en.wikipedia.org/wiki/Continued_fraction_factorization
http://mathworld.wolfram.com/ContinuedFractionFactorizationAlgorithm.html
Notable achievements of this method, first envisaged in 1931, include the factorization of F7, the seventh Fermat number, in 1970 by Morrison and Brillhart.

Saturday 24 November 2007

Rho:x^2-2

Most polynomials work nicely with Pollard Rho.
However, f(x)=x^2 and f(x)=x^2-2 should be avoided.
Here's what John Pollard had to say by way of reason, in 1975 [BIT 15, P.333]:
"(i) that all polynomials x^2+b seem equally good in (1) except that x^2 and x^2-2 should not be used (whatever the starting value x0), the latter for reasons connected with its appearance in the Lucas-Lehmer test for primality of the Mersenne Numbers [3],"
while Knuth has this to say [TAOCP Vol.2 P.386]:
"In those rare cases where failure occurs for large N, we could try using f(x)=x^2+c for some c<>0 or 1. The value c=-2 should also be avoided, since the recurrence x_(m+1) = (x_m)^2-2 has solutions of the form x_m = r^(2^m) + r^-(2^m). Other values of c do not seem to lead to simple relationships mod p, and they should all be satisfactory when used with suitable starting values."

Friday 23 November 2007

2007 chips

I recently came across the following page, which has some interesting speed comparisons of the various types of modern chips...
http://www.tomshardware.com/2007/07/16/cpu_charts_2007/page36.html

Thursday 22 November 2007

Desktop2


This is Lenny's desktop.
It's currently running mprime255OSX on M+2's.
I like the new 'spaces' feature of Leopard, it's something that's been good about Linux desktops for a while.

Wednesday 21 November 2007

LIM (Part 4) - Schonhage-Strassen

The Schonhage-Strassen (SSA) is an asymptotically fast multiplicative algorithm for large integers, developed in 1971.
http://en.wikipedia.org/wiki/Schönhage-Strassen_algorithm
It uses Fast Fourier transforms (FFTs) (more on these at a later date hopefully) and its run-time complexity is of order nlognloglogn. [Note that the FFT must be performed modulo 2^n+1 for a suitable n, but by choosing n large enough this equates to a regular multiplication]
This means that SSA outperforms Karatsuba or Toom-Cook for numbers with tens of thousands of digits or more. An example of its implementation is in GIMPS' Prime95/mprime software. A second example is the recent addition of SSA to the open-source math library GMP.

Tuesday 20 November 2007

LIM (Part 3) - Toom-Cook

Another method of multiplication is called Toom-Cook, first described in 1963.
http://en.wikipedia.org/wiki/Toom-Cook_multiplication
This is basically a generalization of the Karatsuba Method, by splitting the input numbers into multiple parts at a time, rather than just two (as in Karatsuba) or 1 (ie no splitting) (as in classical long multiplication).
Toom-3 (3-way Toom-Cook) reduces 9 multiplications to 5, and runs in order n^(log5/log3) time.

Monday 19 November 2007

WEP-M+2 milestone

The WEP-M+2 Project has just announced that it has reached the milestone of 1000 instances of finding the 12-digit factor of (M+2)2203. Thanks to everyone who has participated so far.
http://bearnol.is-a-geek.com/wanless2/
I estimate (I'm the project admin :) that those 1000 instances equate to about 27 CPU-years (modern CPU cores).

LIM (Part 2) - Karatsuba

A first improvement to long multiplication is the Karatsuba Algorithm.
http://en.wikipedia.org/wiki/Karatsuba_algorithm
This was invented in 1960, and has a time complexity of order n^(log_2(3)).
It relies on the observation that two-digit multiplication can be done with only 3, rather than 4, multiplications classically required. By "dividing and conquering" (ie splitting) the numbers to be multiplied this can be extended to larger numbers.

Sunday 18 November 2007

Large Integer Multiplication (LIM) Part 1

Large Integer Multiplication is an algorithm to multiply two large integers together efficiently.
This is often used by factorization algorithms, both for exponentiation, for example - when it is combined with the Russian Peasant method for extra speed (see earlier) or in fact, often, for division (by multiplying by the inverse of a number as the dividend)
"Long Multiplication" is the obvious, and naive, method, but the time complexity of this is of order n^2 for two n-digit integers. So a number of improvements have been suggested. These will be examined in later parts of this series. For the moment read all about LIM on Wikipedia:
http://en.wikipedia.org/wiki/Multiplication_algorithm

Saturday 17 November 2007

Integer Factorization Records

"Integer factorization records" on Wikipedia has a summary of the current, and recent, state-of-play for the biggest non-trivial numbers yet factored:
http://en.wikipedia.org/wiki/Integer_factorization_records
Note that these have all been achieved with some form of Number Field Sieve, either General (for numbers of no especial form), or Special (for the two Mersenne numbers cited).

Friday 16 November 2007

Msieve v1.30

New version [primarily a bugfix] of msieve (by Jason Papadopoulos) now available.
Download from:
http://www.boo.net/~jasonp/qs.html

The Magic Words are Squeamish Ossifrage


Here is a rather magnificent looking ossifrage [picture from Wikipedia] (though I don't know whether he's especially squeamish! :)
The connection with factorization?
Well the first ever RSA challenge, posed by Martin Gardner in 1977, had this phrase as its encrypted solution - deciphered in 1994 for the $100 prize:
http://en.wikipedia.org/wiki/The_Magic_Words_are_Squeamish_Ossifrage

Thursday 15 November 2007

SQUFOF

Shanks' square forms factorization was devised as an improvement on Fermat's method.

Here is its entry on Wikipedia:
http://en.wikipedia.org/wiki/SQUFOF

and here is its implementation in superfac9:
BigInteger factorizeshanks(BigInteger n) {
BigInteger a = new BigInteger("0");
BigInteger f = new BigInteger("0");
BigInteger h1 = new BigInteger("0");
BigInteger h2 = new BigInteger("0");
BigInteger k = new BigInteger("0");
BigInteger p = new BigInteger("0");
BigInteger pp = new BigInteger("0");
BigInteger q = new BigInteger("0");
BigInteger qq = new BigInteger("0");
BigInteger qqq = new BigInteger("0");
BigInteger r = new BigInteger("0");
BigInteger te = new BigInteger("0");
BigInteger i = new BigInteger("0");
BigInteger count = new BigInteger("0");

k = sqrt(n);

if (fastsquareQ(n)) return k;

a=k; h1=k; h2=ONE; pp=ZERO; qq=ONE; qqq=n; r=ZERO;

for (count=ONE;count.compareTo(TENTHOUSAND)<0;count=count.add(ONE)) {
p=k.subtract(r);
q=qqq.add(a.multiply(pp.subtract(p)));
a=(p.add(k)).divide(q);
r=(p.add(k)).remainder(q);
te=(a.multiply(h1)).add(h2);
h2=h1;
h1=te;
pp=p;
qqq=qq;
qq=q;
te = sqrt(q);
i=i.add(ONE);
if ((i.remainder(TWO).compareTo(ZERO))!=0 || !fastsquareQ(q)) continue;
te=h2.subtract(te);
f=n.gcd(te);
if (f.compareTo(ONE) > 0 && f.compareTo(n) < 0)
return f;
}

return f;
}

Wednesday 14 November 2007

RSA-155

RSA-155 (a 155 digit semiprime) was factored on August 22, 1999 by GNFS.
RSA-155 = 102639592829741105772054196573991675900716567808038066803341933521790711307779
* 106603488380168454820927220360012878679207958575989291522270608237193062808643

Read more about it at Wikipedia:
http://en.wikipedia.org/wiki/RSA-155
and on the official announcement:
http://listserv.nodak.edu/cgi-bin/wa.exe?A2=ind9908&L=nmbrthry&P=1905

Tuesday 13 November 2007

snfspoly

There is also 'snfspoly' - the equivalent of 'phi' for XYYXF composites, which can be difficult to find - search for it on the yahoo XYYXF mailing list, or email me...

Monday 12 November 2007

Phi

Alex Kruppa has written a small (but growing) program, licensed under the GPL, called 'phi', for generating SNFS polynomials for use with GGNFS or msieve. It can currently produce correct polys for cyclotomic numbers (eg Cunningham Project) and 'Homogeneous' Cunninghams. Search for the source code (written in 'C', and using GMP) in the factoring section of the mersenneforum.

Saturday 10 November 2007

What does the term "embarrassingly parallel" mean?

An algorithm that can be run (very) profitably on many threads simultaneously. For example, trial-division, with a different random seed in each thread as the potential factor, to avoid repeating work.

Friday 9 November 2007

FireStream

I imagine these new chips from AMD would be great for embarrassingly parallel apps like most factorization methods...
http://www.pcworld.com/article/id,139413-c,amd/article.html
Here is some more on this series of chips, from Wikipedia:
http://en.wikipedia.org/wiki/AMD_Stream_Processor

Thursday 8 November 2007

Velocity Engine

Apple has produced a sample application, to demonstrate their "Velocity Engine", which uses the vector capability of PPC chips, G4 & G5.
http://en.wikipedia.org/wiki/AltiVec
It happens to be a factorization program! :)
http://developer.apple.com/samplecode/VelEng_Multiprecision/index.html
The program proceeds by trial-division, then rho, and finally ECM - and is quite fast...

Wednesday 7 November 2007

Lenny

Well I've only been and gone and done it! Bought a new MacBook that is, complete with Leopard. It's called "Lenny"
I've already got it running George Woltman's mprime 25.5 for Mac OSX (beta), on medium-sized M+2 numbers, and have high expectations of finding a new factor or two :)

Tuesday 6 November 2007

Richard Brent


Here is a picture of Richard Brent, from his homepage:
http://wwwmaths.anu.edu.au/~brent/
(image linked to in situ)

Monday 5 November 2007

Brent's Method

Brent's improvement (in 1980) to Pollard's rho method is to extend the core idea by generalizing the cycle by powers of two.
Full details here, with due reference to Floyd -
http://web.comlab.ox.ac.uk/oucl/work/richard.brent/pd/rpb051i.pdf

The following from Wikipedia:
Input: n, the integer to be factored; x0, such that 0 ≤ x0 ≤ n; m such that m > 0; and f(x), a pseudo-random function modulo n.
Output: a non-trivial factor of n, or failure.
1. y ← x0, r ← 1, q ← 1.
2. Do:
1. x ← y
2. For i = 1 To r:
1. y ← f(y)
3. k ← 0
4. Do:
1. ys ← y
2. For i = 1 To min(m, r − k):
1. y ← f(y)
2. q ← (q × |x − y|) mod n
3. g ← GCD(q, n)
4. k ← k + m
5. Until (k ≥ r or g > 1)
6. r ← 2r
3. Until g > 1
4. If g = n then
1. Do:
1. ys ← f(ys)
2. g ← GCD(|x − ys|, n)
2. Until g > 1
5. If g = n then return failure, else return g

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function brentFactor(N)
# Initial values x(i) and x(m) for i = 0.
xi := 2
xm := 2
for i from 1 to infinity
# Find x(i) from x(i-1).
xi := (xi ^ 2 + 1) % N
s := gcd(xi - xm, N)
if s <> 1 and s <> N then
return s, N/s
end if
if integralPowerOf2(i) then
xm := xi
end if
end do
end function

Here is its Java implementation from superfac9:

BigInteger factorizebrent(BigInteger n) {
BigInteger k = new BigInteger("1");
BigInteger r = new BigInteger("1");
BigInteger i = new BigInteger("1");
BigInteger m = new BigInteger("1");
BigInteger iter = new BigInteger("1");
BigInteger z = new BigInteger("1");
BigInteger x = new BigInteger("1");
BigInteger y = new BigInteger("1");
BigInteger q = new BigInteger("1");
BigInteger ys = new BigInteger("1");

m=TEN;
r=ONE;
iter=ZERO;
z=ZERO;
y=z;
q=ONE;

do {
x=y;
for (i=ONE;i.compareTo(r)<=0;i=i.add(ONE)) y=((y.multiply(y)).add(THREE)).mod(n);
k=ZERO;
do {
iter=iter.add(ONE);
// System.out.print("iter=" + iter.toString() + '\r');
ys=y;
for (i=ONE;i.compareTo(mr_min(m,r.subtract(k)))<=0;i=i.add(ONE)) {
y=((y.multiply(y)).add(THREE)).mod(n);
q=((y.subtract(x)).multiply(q)).mod(n);
}
z=n.gcd(q);
k=k.add(m);
} while (k.compareTo(r)<0 && z.compareTo(ONE)==0);
r=r.multiply(TWO);
} while (z.compareTo(ONE)==0 && iter.compareTo(TENTHOUSAND)<0);

if (z.compareTo(n)==0) do {
ys=((ys.multiply(ys)).add(THREE)).mod(n);
z=n.gcd(ys.subtract(x));
} while (z.compareTo(ONE)==0);

return z;
}

Achievements of this method include the factorization, in 1980, of the eighth Fermat number:
http://wwwmaths.anu.edu.au/~brent/pub/pub061.html

Sunday 4 November 2007

Pollard P-1 Method

There is also the Pollard "p-1" method, invented in 1974. It relies on Fermat's Little Theorem.
Here is a Wikipedia article about this method:
http://en.wikipedia.org/wiki/Pollard%27s_p_-_1_algorithm

And here is its description on the mersenneforum:
http://mersennewiki.org/index.php/P-1_Factorization_Method

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function pollard_p1(N)
# Initial value 2^(k!) for k = 0.
two_k_fact := 1
for k from 1 to infinity
# Calculate 2^(k!) (mod N) from 2^((k-1)!).
two_k_fact := modPow(two_k_fact, k, N)
rk := gcd(two_k_fact - 1, N)
if rk <> 1 and rk <> N then
return rk, N/rk
end if
end for
end function

Saturday 3 November 2007

A Random Factorization

Below is the factorization of a random 100-digit number.

First in Pari/GP: (which uses Rho, ECM and MPQS)
[NB the time quoted is process-time, rather than real-time]
? #
timer = 1 (on)
? factor(905771525917281232131519213461223147373627632478259763073719184206592688398458994971036043749073482)
time = 12mn, 17,641 ms.
%1 =
[2 1]

[3 1]

[11 1]

[18701 1]

[111977 1]

[122016508135030794072521 1]

[3174449800530489735869567 1]

[16919752823495547077187437987066464785943 1]

...and then with superfac9:
tiggatoo:~/math james$ time java superfac9 < random100d.txt
[905771525917281232131519213461223147373627632478259763073719184206592688398458994971036043749073482]
wanless...
brutep: 2
wanless...
wanless...
brutep: 3
brutep: 11
ecm...
ecm...
aprtcle: 18701
aprtcle: 111977
ecm...
aprtcle: 122016508135030794072521
siqs...
aprtcle: 3174449800530489735869567
aprtcle: 16919752823495547077187437987066464785943
duration: 27490 seconds

Exception in thread "main" java.lang.StringIndexOutOfBoundsException: String index out of range: 0
at java.lang.String.charAt(String.java:558)
at superfac9.main(superfac9.java:192)

real 458m11.888s
user 59m52.469s
sys 0m29.617s

Friday 2 November 2007

Trial Division

Perhaps I should have started with this, but...
Trial division is the simplest and most naive of factorization methods.
It _is_ however, the fastest method for easily eliminating very small factors from composite inputs.
Note also, though, interestingly, that in fact it is also just about the only method capable of finding small (or indeed any) factors of really large numbers of specific algebraic form(s), because those algebraic forms can then be calculated modulo the suspected factor, keeping the required calculations small.

http://en.wikipedia.org/wiki/Trial_division

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function trialDivision(N)
for s from 2 to floor(sqrt(N))
if s divides N then
return s, N/s
end if
end for
end function

Here is its Java implementation (as 'brute' - standing for 'brute force') from superfac9:

BigInteger factorizebrute(BigInteger n) {
BigInteger a = new BigInteger("2");

while (a.compareTo(TENTHOUSAND) < 0 && a.multiply(a).compareTo(n) <= 0) {
if (n.remainder(a).compareTo(ZERO) == 0)
return a;
else
a = a.add(ONE);
}

return n;
}

Thursday 1 November 2007

P+1 p60

This just in... Alex Kruppa is reporting a new, record-size factor, of 60 digits, found with the P+1 method.
Read all about it here:
http://mersenneforum.org/showthread.php?p=117442#post117442

Wave-particle duality

Today I (also) have some musings on the wave-particle duality for you, and its effect on factorization.
Clearly there is a mapping (according to QT), from the wave-to-particle domains, or very equivalently, the time-to-space domains.
As integers get larger, the required precision to express that integer increases - leading to a more efficient representation/manipulation in the wave domain, rather than the particle. Hence FFT-methods for eg Large-integer-multiplication (more on that some other time hopefully).
Note that a so-called quantum computer, would also be operating naturally largely in the wave-domain.
This leads me to the suggestion, that maybe a QC would operate best on native representations of FFTs...

Russian Peasant

'Russian Peasant' is a method for fast exponentiation. It is also known as 'exponentiation by squaring'. The basic idea is to exploit knowledge of the binary expansion of the exponent, by using selective repeated squaring. Depending on the factorization algorithm, and the size of the input number, this can lead to significant performance enhancement.

http://en.wikipedia.org/wiki/Exponentiation_by_squaring

Here is an on-line demo:
http://www.math.umbc.edu/~campbell/NumbThy/Class/BasicNumbThy.html#Modular-PowMod

The classic description is in Knuth (The Art of Computer Programming) 4.6.3
Java code as below:
 static BigInteger power(BigInteger x, BigInteger n, BigInteger mod) { 
// Knuth 4.6.3 - computes x^n modulo mod 
  BigInteger N = new BigInteger("0"); 
  BigInteger Y = new BigInteger("0"); 
  BigInteger Z = new BigInteger("0");
  N=n; 
  Y=one; 
  Z=x;
  while (true) { 
   if (N.remainder(two).compareTo(zero) > 0) { 
    N=N.divide(two); 
    Y=Z.multiply(Y).remainder(mod); 
    if (N.compareTo(zero) == 0) 
     return Y; 
   } 
   else { 
    N=N.divide(two); 
   } 
   Z=Z.multiply(Z).remainder(mod); 
  } 
 }

Wednesday 31 October 2007

History of factorization records by sieving methods


This graph is by Francois Morain. http://algo.inria.fr/seminars/sem00-01/morain.html
(image linked to in situ)
The relevant figure is towards the bottom of the article, in Section 4.

Tuesday 30 October 2007

Penryn

It would appear that Intel's next-generation processors (Penryn) draw up to 40% less power than current CPUs, under full load. :))) Quite welcome I imagine, in the current climate (pun intended :( ).
http://www.tomshardware.com/2007/10/29/intel_penryn_4ghz_with_air_cooling/page14.html
This is no doubt due to moving to the smaller 45nm die resolution. Paul Otellini (CEO) has said that Penryn will ship 12 November 2007...

Monday 29 October 2007

Coincidence?

Is it coincidence that Fermat's method works best with factors close to the square root of the input number, while NFS is particularly good at factoring semiprimes?

Sunday 28 October 2007

XYYXF Project milestone

Congratulations to the XYYXF project, who have just finished their original target of factoring all composites with x and y <=100, after 6 years! The search continues with x and y extended to 150...

http://xyyxf.at.tut.by/news.html#0

Daylight Savings Time

Now the clocks have gone back, all my ECMNet times are correct again ! :)

Msieve v1.29

New version of msieve (by Jason Papadopoulos) now available.
Download from:
http://www.boo.net/~jasonp/qs.html
Announcement at:
http://mersenneforum.org/showthread.php?p=117235#post117235

Fermat's Method

This 400-yr old method is the basis for many modern sieving factorization methods.
It relies on the discovery of an expression for the input number, as the difference between two (integral) squares.

http://en.wikipedia.org/wiki/Fermat%27s_factorization_method

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function fermatFactor(N)
for x from ceil(sqrt(N)) to N
ySquared := x * x - N
if isSquare(ySquared) then
y := sqrt(ySquared)
s := (x - y)
t := (x + y)
if s <> 1 and s <> N then
return s, t
end if
end if
end for
end function

And here is its Java implementation (from superfac9) [with thanks to DT]:

BigInteger factorizefermat(BigInteger n) {
BigInteger a, bSq;
long iterations = 1L;

a = sqrt(n);
if(n.mod(a).compareTo(ZERO) == 0) return a;

while(iterations<10000000) {
bSq = a.multiply(a).subtract(n);
if(fastsquareQ(bSq))
return a.subtract(sqrt(bSq)); // bSq is a square, factorization found
a = a.add(ONE);
iterations++;
}

return n;
}

Saturday 27 October 2007

Mark Manasse


Here is a picture of Mark Manasse (reproduced by kind permission), who is notable for his involvement in the factorizations of F9, RSA-110 and RSA-120 in the early 1990's. In fact he helped develop the Number Field Sieve, which counted these numbers among its early successes.
http://www.std.org/~msm/common/f9paper.pdf
[Incidentally this paper, above, I see, also explains why WE works on prime-powers...basically Fermat's Little Theorem in combination w/ the Binomial Theorem]
http://en.wikipedia.org/wiki/RSA-110
http://en.wikipedia.org/wiki/RSA-120
http://www.std.org/~msm/common/nfspaper.pdf

Friday 26 October 2007

"factorization"

Here's a bit of a Meta-post:
Regarding the spelling/usage of the word "factorization".
First off note that there is a synonym for "factorization", namely "factoring" - however since this also has a completely unrelated meaning (financial of some sort IIRC), I generally tend to use the longer form, unless it would sound particularly ugly.
Similarly the verb can just be "to factor", instead of "to factorize", however in this case the correct (mathematical) meaning is usually clearer, since it will be associated with an object (eg C127_113_36).
In addition (just to add to the confusion :) there is an alternative _spelling_ of "factorization" as "factorisation" (in British English).

Thursday 25 October 2007

Desktop


My desktop is pretty busy at the moment.
Here is a current screen grab.
The windows are (left-to-right):
1) Random-based WEP on (M+2)859433
2) WEP on F7
3) WEP on F8
4) WEP on (M+2)4253
5) msieve (NFS) on C127_113_36 (xyyxf)
6) ECM via ECMNet for xyyxf
7) GGNFS on rsa100
8) GGNFS on C127_113_36
All this on a 1.9GHz G5 (iMac called 'tiggatoo')!

Wednesday 24 October 2007

GGNFS on PPC (update)

Basically, the problem of "not enough polynomials in mpqs", is fixed in the latest CVS code (0.77.1-20060722) available on Sourceforge (with a separate, preliminary, classical sieve stage).
This is the story of how I succeeded in obtaining, compiling and running it on my iMac G5 (with thanks to Mark Rodenkirch, who responded to my query on the yahoo GGNFS mailing list).
Step 1)
Obtain, compile and install a 64-bit GMP archive (PPC 970)
Step 2)
Download the source code of GGNFS via CVS from Sourceforge, thus:
"cvs -d:pserver:anonymous@ggnfs.cvs.sourceforge.net:/cvsroot/ggnfs
checkout -R branch_0"
The source-tree is also available for browsing at:
http://ggnfs.cvs.sourceforge.net/ggnfs/branch_0/src/
Step 3)
Necessary specifically (atm) for PPC [(probably) not for other architectures]:
tweak/hack the code to include its own implementation [supplied, but disabled] of getline() in if.c by,
a) inserting a suitable declaration for the missing getline() function at the top of file if.c
b) remming out the ifdef (and matching endif) of GGNFS_GNU_MISSING lower down, to actually activate the missing implementation of getline()
(I also investigated other methods of achieving same eg trying to define GGNFS_GNU_MISSING elsewhere to achieve the same effect, but this was the first method that worked.)
Step 4)
Compile with "make ppc_970"
Step 5)
Note that the default implementation of factLat.pl, by which factorizations are actually run, already has the correct path to the binaries, so does not need adjusting, unlike some previous versions.

Tuesday 23 October 2007

Fermat number factor search

For the very latest on the search for Fermat number factors, checkout George Woltman's page at:
http://www.mersenne.org/ecmf.htm
Prime95/mprime includes an ECM factoring ability, easily controlled by a single ini file, and is particularly well-suited (and fast!) to looking for factors of larger Fermat (and Mersenneplustwo ;-) numbers. George Woltman even supplies a file of already-known factors of Fermat numbers, that mprime can read, to avoid the time spent re-discovering these (or you can start from scratch without reading-in that file, as a double-check that your set-up is finding factors successfully - this won't affect your ability to also find new factors)

Monday 22 October 2007

Euclidean Algorithm

The Euclidean Algorithm is vital to many factorization algorithms.
It is a very fast method for finding the gcd (greatest common divisor) of two numbers, without actually performing the factorization of either.
This page has pretty much all you need to know about it:
http://en.wikipedia.org/wiki/Euclidean_algorithm

Sunday 21 October 2007

"Monte Carlo" method(s)

What is a "Monte Carlo" (factorization) method?
Essentially one where luck can play a major part - hence the rather exciting-sounding connection with gambling! This is the case whenever the process is seeded by some (pseudo-or-otherwise) random number generation (eg for the benefit of parallelization), or when there is inherent unpredictability in the algorithm itself.
Note also, that for the Rho method, for example, this Monte Carlo nature of the algorithm results in performance gains caused by a statistical effect akin to the so-called "Birthday Paradox". This latter, in turn, is the term given to the surprisingly high observed probability of at least one coincidence of two different random variables in the same range (eg two or more people's birthdays)
The following Wikipedia link talks about Monte Carlo effect in general:
http://en.wikipedia.org/wiki/Monte_Carlo_method
and here is a link to the Birthday Paradox:
http://en.wikipedia.org/wiki/Birthday_Paradox

Saturday 20 October 2007

Rho method

Today I am going to be talking a little bit about Pollard's Rho algorithm, which, for such an apparently simple algorithm, is surprisingly good at finding moderate-size factors of general numbers.

Some reference links:
http://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
[as I mentioned briefly once before]
http://mersennewiki.org/index.php/Rho_Factorization_Method

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function pollardRho(N)
# Initial values x(i) and x(2*i) for i = 0.
xi := 2
x2i := 2
do
# Find x(i+1) and x(2*(i+1))
xiPrime := xi ^ 2 + 1
x2iPrime := (x2i ^ 2 + 1) ^ 2 + 1
# Increment i: change our running values for x(i), x(2*i).
xi := xiPrime % N
x2i := x2iPrime % N
s := gcd(xi - x2i, N)
if s <> 1 and s <> N then
return s, N/s
end if
end do
end function


The Rho method is also the main method used by GMP's demo factorization sample program (or was last time I looked) - immediately below:

void
factor_using_pollard_rho (mpz_t n, int a_int, unsigned long p)
{
mpz_t x, x1, y, P;
mpz_t a;
mpz_t g;
mpz_t t1, t2;
int k, l, c, i;

if (flag_verbose)
{
printf ("[pollard-rho (%d)] ", a_int);
fflush (stdout);
}

mpz_init (g);
mpz_init (t1);
mpz_init (t2);

mpz_init_set_si (a, a_int);
mpz_init_set_si (y, 2);
mpz_init_set_si (x, 2);
mpz_init_set_si (x1, 2);
k = 1;
l = 1;
mpz_init_set_ui (P, 1);
c = 0;

while (mpz_cmp_ui (n, 1) != 0)
{
S2:
if (p != 0)
{
mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
}
else
{
mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
}
mpz_sub (t1, x1, x); mpz_mul (t2, P, t1); mpz_mod (P, t2, n);
c++;
if (c == 20)
{
c = 0;
mpz_gcd (g, P, n);
if (mpz_cmp_ui (g, 1) != 0)
goto S4;
mpz_set (y, x);
}
S3:
k--;
if (k > 0)
goto S2;

mpz_gcd (g, P, n);
if (mpz_cmp_ui (g, 1) != 0)
goto S4;

mpz_set (x1, x);
k = l;
l = 2 * l;
for (i = 0; i < k; i++)
{
if (p != 0)
{
mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
}
else
{
mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
}
}
mpz_set (y, x);
c = 0;
goto S2;
S4:
do
{
if (p != 0)
{
mpz_powm_ui (y, y, p, n); mpz_add (y, y, a);
}
else
{
mpz_mul (y, y, y); mpz_add (y, y, a); mpz_mod (y, y, n);
}
mpz_sub (t1, x1, y); mpz_gcd (g, t1, n);
}
while (mpz_cmp_ui (g, 1) == 0);

if (!mpz_probab_prime_p (g, 3))
{
do
{
mp_limb_t a_limb;
mpn_random (&a_limb, (mp_size_t) 1);
a_int = (int) a_limb;
}
while (a_int == -2 || a_int == 0);

if (flag_verbose)
{
printf ("[composite factor--restarting pollard-rho] ");
fflush (stdout);
}
factor_using_pollard_rho (g, a_int, p);
break;
}
else
{
mpz_out_str (stdout, 10, g);
fflush (stdout);
fputc (' ', stdout);
}
mpz_div (n, n, g);
mpz_mod (x, x, n);
mpz_mod (x1, x1, n);
mpz_mod (y, y, n);
if (mpz_probab_prime_p (n, 3))
{
mpz_out_str (stdout, 10, n);
fflush (stdout);
fputc (' ', stdout);
break;
}
}

mpz_clear (g);
mpz_clear (P);
mpz_clear (t2);
mpz_clear (t1);
mpz_clear (a);
mpz_clear (x1);
mpz_clear (x);
mpz_clear (y);
}


This Java code is taken directly from superfac9:

BigInteger factorizerho(BigInteger n) {
BigInteger loop = new BigInteger("1");
BigInteger x = new BigInteger("5");
BigInteger y = new BigInteger("2");

while (n.gcd(x.subtract(y)).compareTo(ONE) == 0 && loop.compareTo(TENTHOUSAND) < 0) {
x = x.multiply(x).add(ONE).mod(n);
x = x.multiply(x).add(ONE).mod(n);
y = y.multiply(y).add(ONE).mod(n);
loop = loop.add(ONE);
}

return n.gcd(x.subtract(y));

}


Note that there are several variants (given different start-point values for x and y for example), but the basic "Monte-Carlo" principle remains the same.

Notable successes with the Rho method or its variants include Brent and Pollard's factorization of the eighth Fermat number in 1980. The 16-digit factor of F8 took about 2 hours to find on a UNIVAC 1100/42, apparently.

Friday 19 October 2007

F12-WEP

And here is the WEP algorithm running (unsuccessfully) for an hour or so looking for new factors of F12:

tiggatoo:~/math/wec james$ time ./factorize3.gmp -P4096 -T1000 114689 26017793 63766529
P4096 T1000

real 59m16.287s
user 12m32.262s
sys 0m4.508s

Thursday 18 October 2007

M521^2 by WE

Here is a factorization illustrating the automatic algebraic factoring capability of superfac9, due to WE algorithm:

tiggatoo:~/math james$ java superfac9
f(2^521-1)^2

[47125446914534694131579097993419809976955095716785201420286055195012674566357244479460731079205201122720511132925006540350105785156086431086764996857554304847155991333706718342307167456986269662311038377104760933477381254100896222805785374204495333936040246318307567782851014765052850751581472024524956029996236801]
wanless...
aprtcle: 6864797660130609714981900799081393217269435300143305409394463459185543183397656052122559640661454554977296311391480858037121987999716643812574028291115057151
aprtcle: 6864797660130609714981900799081393217269435300143305409394463459185543183397656052122559640661454554977296311391480858037121987999716643812574028291115057151
duration: 115 seconds

Wednesday 17 October 2007

Near-repdigit factorization(s)

And speaking of repunits (or repdigits - the factorization of these latter reduce to the former), there is also Makoto Kamada's page of factorizations of near-repdigits.
http://homepage2.nifty.com/m_kamada/math/factorizations.htm
Clearly this enlarges the options for factoring candidates significantly.

Tuesday 16 October 2007

Repunit factorization

Here is a link to Yousuke Koide's page on repunit factorization(s). Repunits are numbers specifically whose decimal representation consists of all 1's.
http://www.h4.dion.ne.jp/~rep/
It would appear that the lowest repunit whose prime factorization is yet unknown, is R239 (see also Torbjorn Granlund's page at Swox).
http://swox.com/~tege/fac10m.txt
[Note that Swox is the organization behind the widely-used GMP math library
http://gmplib.org/]

Monday 15 October 2007

GGNFS bounce!


GGNFS is cool! (and deservedly a standard reference implementation). But it would appear I have been a victim of the GGNFS "bounce", as Bob Backstrom would call it! I had almost collected enough relations to proceed to the linear algebra stage of my factorization attempt of C127_113_36 (for the xyyxf project), when this happened. Apparently, discarding heavy relations speeds up the later stage, but it's going to take me a while to replace those discarded relations. Bob assures me that this bounce will not happen again, and also suggests (as in fact I wondered) that it may well be possible to proceed to the LA immediately with msieve. However, since I'm fairly new at this (not having run even straight GGNFS before on a number this size), for now I think I've decided to wait, watch, and hopefully learn...

Sunday 14 October 2007

Benchmarks: Java vs C (Addendum)

Of course Assembler can be even faster than C. George Woltman's big number library (gwnum) - part of Prime95/mprime - for example, seems very fast, perhaps at least as fast again as GMP, particularly for extrabig numbers.
Gwnum also uses advanced FFT techniques, and because it is partly written in assembler, not surprisingly, it is not portable - only being available for the most common platform, x86.
It would be interesting to rewrite my WEP application using gwnum...
http://www.mersenne.org/freesoft.htm
http://en.wikipedia.org/wiki/Fft

Saturday 13 October 2007

Benchmarks: Java vs C (Part 4)

Here are the promised benchmarks, running random-based WEP on input numbers of, respectively, 100, 1000 and 10000 bits.
If anyone would care to (and has the patience to) extend the tables, feel free to send me results!
Result: Java is, very approximately, an order of magnitude, only, slower than C.

tiggatoo:~/math/wec james$ time java we2tpr2 < P109.in.2.txt
[109]
[1]
base#0
elapsed=0s factor=3 A=645264517869175453037177387327
base#10746
elapsed=43s factor=104124649 A=245093880973395622884386051056621
2077756847362348863128179
duration: 43 seconds


real 0m46.195s
user 0m5.185s
sys 0m0.510s

tiggatoo:~/math/wec james$ time ./factorize3.gmp -P109 -T10748 3
P109 T10748

real 0m7.176s
user 0m0.793s
sys 0m0.011s

tiggatoo:~/math/wec james$ time java we2tpr2 < P1279.in.2.txt
[1279]
[1]
base#0
elapsed=2s factor=3 A=530978962560898025744780713800683
base#11
elapsed=13s factor=706009 A=1243358263022091534365334367195221
^Cse#100
real 2m34.366s
user 0m17.542s
sys 0m0.222s

tiggatoo:~/math/wec james$ time ./factorize3.gmp -P1279 -T114 3 706009
P1279 T114

real 0m29.889s
user 0m3.344s
sys 0m0.022s

tiggatoo:~/math/wec james$ time java we2tpr2 < P11213.in.2.txt
[11213]
[1]
base#0
elapsed=1510s factor=3 A=93613194278166973511368466002199
^Cse#2
real 45m19.926s
user 5m14.923s
sys 0m2.640s

tiggatoo:~/math/wec james$ time ./factorize3.gmp -P11213 -T3 3
P11213 T3

real 5m44.733s
user 0m40.864s
sys 0m0.247s

Friday 12 October 2007

Benchmarks: Java vs C (Part 3)

Here is the (nearly) same program (factorize3.mac.c) coded in C (also copyright JGW). It needs to be compiled with -lgmp, the GMP math library.

/* Factoring with WEP method using random base(s).
*/

#include
#include
#include

#include

#include "gmp.h"

int flag_verbose = 0;

int
factor_using_random_wep (mpz_t s, unsigned long p, unsigned long T)
{
mpz_t V, v1, v2, v3, b, bb, bbb, bbbb, wanless;
long numtrials;
int flag;

gmp_randstate_t rstate;



if (flag_verbose)
{
printf ("[wep ");
printf ("s=");
mpz_out_str (stdout, 10, s);
printf ("\tp=");
printf ("%ld", p);
printf ("\tT=");
printf ("%ld", T);
printf ("]\n\r");
fflush (stdout);
}


if (mpz_probab_prime_p (s, 3))
{
printf ("P%ld\tT0\t", p);
mpz_out_str (stdout, 10, s);
fflush (stdout);
flag=2;
return(flag);
}



gmp_randinit (rstate, GMP_RAND_ALG_LC, 128);

{
#if HAVE_GETTIMEOFDAY
struct timeval tv;
gettimeofday (&tv, NULL);
gmp_randseed_ui (rstate, tv.tv_sec + tv.tv_usec);
#else
time_t t;
time (&t);
gmp_randseed_ui (rstate, t);
#endif
}






numtrials=0;

mpz_init_set_si (b, 1);
mpz_init_set_si (bb, 1);
mpz_init_set_si (bbb, 1);
mpz_init_set_si (bbbb, 1);
mpz_init_set_si (V, 1);
mpz_init_set_si (v1, 1);
mpz_init_set_si (v2, 1);
mpz_init_set_si (v3, 1);
mpz_init_set_si (wanless, 2);

while (mpz_cmp (wanless, s) < 0)
mpz_mul_ui (wanless, wanless, 2);


while (numtrials < T && (mpz_cmp_ui (V, 1) == 0 || mpz_cmp (V, s) == 0))
{


mpz_urandomb (b, rstate, 100L);
mpz_mul_ui(bb, b, p);
mpz_mul_ui(bbb, bb, 2);
mpz_add_ui(bbbb, bbb, 1);


mpz_powm (v1, bbbb, wanless, s);
mpz_powm (v2, v1, s, s);
mpz_sub (v3, v2, v1);
mpz_gcd (V, s, v3);

if (flag_verbose)
if (numtrials%1000 == 0)
{
printf ("numtrials=%ld\tb=", numtrials);
mpz_out_str (stdout, 10, b);
printf ("\r");
fflush (stdout);
}

numtrials++;

}

if (flag_verbose)
printf("\n\r");


if (mpz_cmp_ui (V, 1) > 0 && mpz_cmp (V, s) < 0)
{
printf ("P%ld\tT%ld\t", p, numtrials);
mpz_out_str (stdout, 10, V);
printf("\tbase=");
mpz_out_str (stdout, 10, bbbb);
printf ("\n");
fflush (stdout);

flag=3;
}

else
{
printf ("P%ld\tT%ld\n", p, T);
flag=0;
}



mpz_clear (b);
mpz_clear (bb);
mpz_clear (bbb);
mpz_clear (bbbb);
mpz_clear (V);
mpz_clear (v1);
mpz_clear (v2);
mpz_clear (v3);
mpz_clear (wanless);

return (flag);

}

main (int argc, char *argv[])
{
mpz_t r, s, t, F, f;
unsigned long p, T;
int i;
int flag;

if (argc > 1 && !strcmp (argv[1], "-v"))
{
flag_verbose = 1;
argv++;
argc--;
}

mpz_init (r);
mpz_init (s);
mpz_init (t);
mpz_init (F);
mpz_init (f);

mpz_set_ui (F, 1);
mpz_set_ui (f, 1);

if (argc > 1)
{
p = 0;
for (i = 1; i < argc; i++)
{
if (!strncmp (argv[i], "-M", 2))
{
p = atoi (argv[i] + 2);
mpz_set_ui (t, 1);
mpz_mul_2exp (t, t, p);
mpz_sub_ui (t, t, 1);
}
else if (!strncmp (argv[i], "-P", 2))
{
p = atoi (argv[i] + 2);
mpz_set_ui (t, 1);
mpz_mul_2exp (t, t, p);
mpz_add_ui (t, t, 1);
}
else if (!strncmp (argv[i], "-T", 2))
{
T = atoi (argv[i] + 2);
}
else
{
mpz_set_str (f, argv[i], 0);
mpz_mul (F, F, f);
}
}


mpz_mod (r, t, F);

if (mpz_cmp_si (r, 0) != 0) {
printf ("Wrong known factors!\n");
flag =1;
}
else {

mpz_div (s, t, F);

flag = factor_using_random_wep (s, p, T);

}
}

mpz_clear(r);
mpz_clear(s);
mpz_clear(t);
mpz_clear(F);
mpz_clear(f);

exit (flag);
}

Thursday 11 October 2007

Benchmarks: Java vs C (Part 2)

Here is the program, we2tpr2.java, coded in Java (copyright JGW, all rights reserved :) ). Note how elegant, and convenient, the BigInteger facility is - this has been a standard feature of Java from its inception. Java is almost unique in that respect - most other languages rely on add-in libraries to provide big number capability.

import java.math.BigInteger;
import java.util.Random;
import java.util.Date;

public class we2tpr2 {
static BigInteger zero = new BigInteger("0");
static BigInteger one = new BigInteger("1");
static BigInteger two = new BigInteger("2");
static BigInteger hundred = new BigInteger("100");
static BigInteger thousand = new BigInteger("1000");

static BigInteger n = new BigInteger("0");
static BigInteger p = new BigInteger("0");
static BigInteger known = new BigInteger("0");
static BigInteger numtrials = new BigInteger("0");
static String s;

static int olds1 = 0;
static int s1 = 0;

static Date d;
static long starttime;
static long finishtime;
static long duration;




public static void main (String args[])
throws java.io.IOException {

we2tpr2 we2tpr2inst = new we2tpr2();

char c;
String sInput;

StringBuffer sbInput = new StringBuffer("");

while ((c = (char)System.in.read()) != '\n' && c != '\r')
sbInput.append(c);
System.in.read();
sInput = sbInput.toString().trim();


if (sInput.charAt(0) == 'f' || sInput.charAt(0) == 'F') {
s = sInput.substring(1).trim();

s1 = 0;
olds1 = 0;
p = we2tpr2inst.eval(s);

System.out.println('[' + p.toString() + ']');
}
else {
p = new BigInteger(sInput);
}



sbInput = new StringBuffer("");

while ((c = (char)System.in.read()) != '\n' && c != '\r')
sbInput.append(c);
System.in.read();
sInput = sbInput.toString().trim();


if (sInput.charAt(0) == 'f' || sInput.charAt(0) == 'F') {
s = sInput.substring(1).trim();

s1 = 0;
olds1 = 0;
known = we2tpr2inst.eval(s);

System.out.println('[' + known.toString() + ']');
}
else {
known = new BigInteger(sInput);
}



n = we2tpr2inst.mersenneplustwo(p);
if (n.remainder(known).compareTo(zero) > 0) {
System.out.println("Wrong known factors!");
return;
}
else
n = n.divide(known);

d = new Date();
starttime = d.getTime();

we2tpr2inst.factorize(n, p);

d = new Date();
finishtime = d.getTime();

duration = (finishtime-starttime)/1000;

System.out.println("duration: " + duration + " seconds");

System.out.println();
}


public BigInteger mersenneplustwo(BigInteger p) {
BigInteger i = new BigInteger("0");

n = two;
for (i = one; i.compareTo(p) < 0; i = i.add(one))
n = n.multiply(two);

n = n.add(one);

return n;
}


public boolean factorize(BigInteger n, BigInteger p) {
boolean prime = false;
BigInteger numtested = new BigInteger("0");
BigInteger T = new BigInteger("1");
BigInteger b = new BigInteger("1");
BigInteger A = new BigInteger("2");
BigInteger wanless = new BigInteger("2");

if (n.isProbablePrime(1000)) {
prime = true;
System.out.println(n);
return prime;
}

// workaround - apparent java bug in modPow - JGW
if (n.compareTo(two) < 0)
return false;

if (n.remainder(two).compareTo(zero) == 0) {
System.out.println(two.toString());
return true;
}
// end workaround

while (wanless.compareTo(n) < 0)
wanless = wanless.multiply(two);

Random r = new Random();

numtested = zero;

while (T.compareTo(one) == 0 || T.compareTo(n) == 0) {
// changed JW 2005-3-23
A = new BigInteger(hundred.intValue(), r);
A = (A.multiply(two).multiply(p)).add(one);

// added JGW 2006-06-09
System.out.print("base#" + numtested + '\r');



// changed DT 2005-2-20
b = A.modPow(wanless, n);
T = n.gcd(b.modPow(n, n).subtract(b));

numtested = numtested.add(one);
}





if (T.compareTo(one) > 0 && T.compareTo(n) < 0) {

d = new Date();
finishtime = d.getTime();

duration = (finishtime-starttime)/1000;


System.out.println();
System.out.println("elapsed=" + duration + "s" + '\t' + "factor=" + T.toString() + '\t' + "A=" + A.toString() + '\t');
factorize(n.divide(T), p);
}

return prime;

}


public BigInteger evalRand(char op, BigInteger oldn) {
BigInteger n = new BigInteger("1");


switch (op) {
case 'r':
case 'R':
Random r = new Random();

n = new BigInteger(oldn.intValue(), r);
break;

default:
n = oldn;
break;
}
return n;
}

public BigInteger evalFact(BigInteger oldn, char op) {
BigInteger n = new BigInteger("1");
BigInteger i = new BigInteger("1");
BigInteger j = new BigInteger("1");
boolean prime = true;


switch (op) {
case '!':
for (i = one; i.compareTo(oldn) <= 0; i = i.add(one))
n = n.multiply(i);
break;

case '#':
for (i = one; i.compareTo(oldn) <= 0; i = i.add(one)) {
prime = true;
for (j = two; (prime == true) && (j.multiply(j).compareTo(i) <= 0); j = j.add(one))
if (i.remainder(j).compareTo(zero) == 0)
prime = false;
if (prime == true)
n = n.multiply(i);
}
break;

default:
n = oldn;
break;
}
return n;
}

public BigInteger evalPower(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '^':
n = oldn.pow(n1.intValue());
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger evalProduct(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '*':
n = oldn.multiply(n1);
break;
case '/':
n = oldn.divide(n1);
break;
case '%':
n = oldn.remainder(n1);
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger evalSum(BigInteger oldn, BigInteger n1, char op) {
BigInteger n = new BigInteger("0");

switch (op) {
case '+':
n = oldn.add(n1);
break;
case '-':
n = oldn.subtract(n1);
break;
default:
n = n1;
break;
}
return n;
}

public BigInteger eval(String s) {
BigInteger oldn0 = new BigInteger("0");
BigInteger oldn1 = new BigInteger("0");
BigInteger oldn2 = new BigInteger("0");
BigInteger n = new BigInteger("0");

char oldop0 = 0;
char oldop1 = 0;
char oldop2 = 0;
char op = 0;

while (s1 < s.length()) {
switch (s.charAt(s1)) {
case '(':
case '[':
case '{':
s1++;
n = eval(s);
break;

case '0':
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
n = readNum(s);
break;

default:
break;
}
if (s1 < s.length()) {
switch (s.charAt(s1)) {

case ')':
case ']':
case '}':
case '!':
case '#':
case 'r':
case 'R':
case '^':
case '*':
case '/':
case '%':
case '+':
case '-':
op = s.charAt(s1);
s1++;
break;

default:
break;
}
}
else
op = 0;


switch (op) {
case 0:
case ')':
case ']':
case '}':
n = evalPower(oldn2, n, oldop2);
n = evalProduct(oldn1, n, oldop1);
n = evalSum(oldn0, n, oldop0);
return n;

case '!':
case '#':
n = evalFact(n, op);
break;

case 'r':
case 'R':
n = readNum(s);
n = evalRand(op, n);
break;

case '^':
n = evalPower(oldn2, n, oldop2);
oldn2 = n;
oldop2 = op;
break;

case '*':
case '/':
case '%':
n = evalPower(oldn2, n, oldop2);
oldop2 = 0;
n = evalProduct(oldn1, n, oldop1);
oldn1 = n;
oldop1 = op;
break;

case '+':
case '-':
n = evalPower(oldn2, n, oldop2);
oldop2 = 0;
n = evalProduct(oldn1, n, oldop1);
oldop1 = 0;
n = evalSum(oldn0, n, oldop0);
oldn0 = n;
oldop0 = op;
break;
default:
break;
}
}
return n;
}

public BigInteger readNum(String s) {
BigInteger n = new BigInteger("0");

olds1 = s1;
while (s1 < s.length() && Character.isDigit(s.charAt(s1)))
s1++;
n = new BigInteger(s.substring(olds1, s1));

return n;
}


}

Wednesday 10 October 2007

Benchmarks: Java vs C (Part 1)

Over the next few days I am going to be benchmarking a small factorization program (random-based WEP algorithm) written (approximately) identically in Java (using Java's built-in BigInteger feature) and in C (using the GMP math library), and seeing how the two compare for various size inputs.
But first some background information about Java:
http://en.wikipedia.org/wiki/Java_%28programming_language%29
and a picture (from Wikipedia) of the Java mascot, Duke:

and (also from Wikipedia) of the creator of Java (at Sun Microsystems), James Gosling:

Tuesday 9 October 2007

Fermat number factorizers

Here is a nice page detailing historical (and recent) figures in the search for factors of Fermat numbers. Note that the first computer calculations are dated around 1953.
http://www.fermatsearch.org/history.php

Monday 8 October 2007

What is a "trapdoor" function?

Term first coined by Diffie and Hellman in 1976, refers to a function that is easy to compute forwards, yet very difficult to invert, ie compute backwards. Large integer multiplication/factorization serves this purpose for RSA.
see link:
http://en.wikipedia.org/wiki/Trapdoor_function

Sunday 7 October 2007

Cryptography

I am currently listening (for free, on iTunes - search on "math 55") to Professor James Demmel's Discrete Mathematics series of lectures at University of California, Berkeley - highly recommended (so far). Anyway, Lecture 18 audio is relevant to today's post, which is on public-key cryptography, and, specifically RSA. The reason cryptography is involved is because RSA relies on the difficulty of factoring large semiprimes for keeping encoded messages secret. Some links:
http://en.wikipedia.org/wiki/Public-key_cryptography
http://en.wikipedia.org/wiki/Rsa
http://www.cs.berkeley.edu/~demmel/ma55_Fall07/LectureNotes/Lecture_15_Oct_03.txt