—#!/usr/bin/perl -w
die
'usage: ./shor.pl [number to factor]'
unless
@ARGV
;
use
strict;
use
warnings;
$Quantum::Entanglement::destroy
= 0;
my
$num
=
$ARGV
[0];
# do some early die'ing
die
"$num is a multiple of two, here I am, brain the size..."
unless
$num
%2;
die
"$num is a non-integer, I only have whole numbers of fingers"
unless
$num
==
int
(
$num
);
die
"$num is less than 15"
unless
$num
>= 15;
"Performing initial classical steps:\n"
;
# work out q value
my
$q_power
=
int
(2*
log
(
$num
) /
log
(2)) +1;
my
$q
= 2 **
$q_power
;
# pick some x so that x is coprime to n.
my
$x
;
do
{
$x
=
int
(
rand
$num
) + 1;
}
until
(
$num
%
$x
!= 0 and
$x
> 2);
#ok.. so this misses the point slightly
"Using q:$q, x:$x\nStarting quantum steps\n"
;
# fill up a register with integers from 1..q
my
$prob
= 1/
sqrt
(
$q
);
my
$register1
= entangle(
map
{
$prob
,
$_
} (0..
$q
-1));
# apply transformation F = x**|a> mod n, store in register 2
# (need to do a p_func to avoid overflow while **)
sub
power_mod {
my
(
$state
,
$x1
,
$num1
) =
@_
;
my
$rt
= 1;
return
1
if
$state
== 0;
return
1
if
$state
== 1;
for
(1..
$state
) {
$rt
= (
$rt
*
$x1
) %
$num1
;
}
return
$rt
;
}
"Performing F = x**|a> mod n\n"
;
my
$register2
= p_func(\
&power_mod
,
$register1
,
$x
,
$num
);
# We now observe $register2, thus partially collapsing reg1
my
$k
=
"$register2"
;
"\$register2 collapsed to $k\n"
;
"Finding period of F (this is where you wish for a QCD)\n"
;
# take a ft of the amplitudes of reg1, placing result in reg3
my
$register3
= QFT(
$register1
);
my
$lqonr
=
"$register3"
;
# observe, this must be multiple of q/r
if
(
$lqonr
== 0) {
"Got period of '0', halting\n"
;
exit
(0);
}
my
$period
=
int
(
$q
/
$lqonr
+ 0.5);
# rounding
"Period of F = x**|a> mod n is $period\n"
;
# now given the period, we need to work out the factor of n
# work out the two thingies:
if
(
$period
% 2 != 0) {
"$period is not an even number, doubling to"
;
$period
*=2;
" $period\n"
;
}
my
$one
=
$x
**(
$period
/2) -1;
my
$two
=
$x
**(
$period
/2) +1;
# one and two must have a gcd in common with n, which we now find...
"$one * $two and $num might share a gcd (classical step)\n"
;
my
(
$max1
,
$max2
) = (1,1);
for
(2..
$num
) {
last
if
$_
>
$num
;
unless
((
$num
%
$_
) || (
$one
%
$_
)) {
$max1
=
$_
;
}
unless
((
$num
%
$_
) || (
$two
%
$_
)) {
$max2
=
$_
;
}
}
"$max1, $max2 could be factors of $num\n"
;
__END__;
=head1 NAME
shor - A short demonstration of Quantum::Entanglement
=head1 SYNOPSIS
./shor.pl [number to factor (>14)]
=head1 DESCRIPTION
This program implements Shor's famous algorithm for factoring numbers. A
brief overview of the algorithm is given below.
=head2 The important maths
Given a number B<n> which we are trying to factor, and some other number
which we have guessed, B<x>, we can say that:
x**0 % n == 1 (as x**0 = 1, 1 % n =1)
There will also be some other number, B<r> such that
x**r % n == 1
or, more specifically,
x**(kr) % n ==1
in other words, the function
F(a) = x**a % n
is periodic with period B<r>.
Now, starting from
x**r = 1 % n
x**(2*r/2) = 1 % n
(x**(r/2))**2 - 1 = 0 % n
and, if r is an even number,
(x**(r/2) - 1)*(x**(r/2) + 1) = 0 mod n
or in nice short words, the term on the left is an integer multiple of B<n>.
So long as x**(r/2) != +-1, at least one of the two brackets on the left
must share a factor with B<n>.
Shor's alorithm provides a way to find the periodicity of the function F
and thus a way to calculate two numbers which share a factor with n, it
is then easy to use a classical computer to find the GCD and thus a
factor of B<n>.
=head1 The steps of the algorithm
=head2 1. Remove early trivial cases
We have efficient classical methods for finding that 2 is a factor of 26,
so we do not need to use this method for this.
=head2 2. Pick an integer
Chose a number B<q> so that C<n**2 <= q <= 2n**2>, this is done on a
classical computer. (This is the size we will use for our quantum register.)
=head2 3. Select at random a number coprime to n
Think of some number less than B<n> so that B<n> and B<x> do not share
a common factor (if they do, we already know the answer...).
=head2 4. Fill a quantum register with integers from 0..q-1
This is where we create our first entangled variable, and is the first
non-classical step in this algorithm.
=head2 5. Calculate F, store in a second register
We now calculate C< F(a) = x**a % n> where a represents the superposition
of states in our first register, we store the result of this in our
second register.
=head2 6. Look at register2
We now look at the value of register two and get some value B<k>, this forces
register1 into
a state which can only collapse into values satisfying the equation
x**a % n = k
The probability amplitudes for the remaining states are now all equal to zero,
note that we have not yet looked directly at register1.
=head2 7. Find period of register1
We now apply a fourier transform to the amplitudes of the states in
register1, storing the result as the probability amplitudes for a new
state with the values of register1. This causes there to be a high
probability that the register will collapse to a value which is some
multiple of C<q/r>.
=head2 8. Observe register1
We now observe register1, and use the result to calculate a likely value
for B<r>. From this we can easily calculate two numbers, one of which
will have a factor in common with n, by applying an efficient classical
algoirthm for finding the greatest common denominator, we will be able
to find a value which could be a factor of B<n>.
=head1 Things to remember
This algorithm does not claim to produce a factor of our number the first
time that it is run, there are various conditions which will cause it
to halt mid-way, for instance, the FT step can give a result of 0 which
is clearly useless. The algorithm is better than any known classical one
because the expectation value of the time required to get a correct answer
is still O(n).
This also cannot factor a number which is prime (it being, as it were, prime)
and also cannot factor something which is a prime power (25, say).
=head1 COPYRIGHT
This code is copyright (c) Alex Gough (alex@rcon.org )2001. This is
free software, you may use, modify and redistribute it under the same
terms as Perl itself.
=head1 BUGS
This is slow, being run on classical computers, ah well.
=cut