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;
}


}

No comments: