CS 596 Quantum Computing


The Feynman Computer Simulator

Main Program Ket Vector Generator Qbit Evaluation

RSA Algorithm Example

RSA Example s^p mod n function

Shor's Algorithm Simulation

Shor's Algorithm Continued Fraction Expansion

 


Main Progam


%main program

%n=#of qbits

n=4

I=sqrt(-1);

BV=Basis(n);

%input the Hamiltonian

%A=[0 1;1 0]^(1/2);

A=1/sqrt(2)*[1 1;1 -1];

%H1=kron(eye(2^(n-1)),A);

H1=kron(kron(eye(4),A),A);

%Feynman's Hamiltonian

c=[0 0;1 0]; a=[0 1;0 0]; E=eye(2);

c1=kron(kron(kron(E,c),E),E);

a0=kron(kron(kron(a,E),E),E);

c2=kron(kron(kron(E,E),c),E);

a1=kron(kron(kron(E,a),E),E);

%H=c1*a0*H1+c2*a1*H1;

H=c1*a0*H1;

H=H+H';

%calculate the Unitary matrix

deltat=1;

U=expm(-I*deltat*H);

%U=U^(1/2);

%Initial Vector

p=zeros(2^n,1);

p(9)=1;

%Evolve the system

for i=1:5

Up=U*p;

T=[real(Up),imag(Up), BV, Up.*conj(Up)];

fprintf(1,' time = %3.0f\n',i)

fprintf(1,'%10.4f + %7.4f I ket[%1.0f %1.0f %1.0f %1.0f] %6.4f\n',T.')

%read the cursor

Up=readq(1,Up);

Up=readq(2,Up);

%Up=readq(3,Up);

%Up=readq(4,Up);

T=[real(Up),imag(Up), BV, Up.*conj(Up)];

fprintf(1,' Read the Cursor time = %3.0f\n',i)

fprintf(1,'%10.4f + %7.4f I ket[%1.0f %1.0f %1.0f %1.0f] %6.4f\n',T.')

p=Up;

end


Top

Ket Vector Generator


function BV=Basis(n)

%generate the basis vectors

m=2^n;

BV=zeros(m,n); P=[0;1];

BV(:,1)=kron(P,ones(2^(n-1),1));

for i=2:(n-1)

BV(:,i)=kron(kron(ones(2^(i-1),1),P),ones(2^(n-i),1));

end

BV(:,n)=kron(ones(2^(n-1),1),P);


Top

Qbit Evaluation


function W=readq(i,Up)

A=[1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0

1 1 1 1 0 0 0 0 1 1 1 1 0 0 0 0

1 1 0 0 1 1 0 0 1 1 0 0 1 1 0 0

1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0];

M=[A(i,:);ones(1,16)-A(i,:)];

pr=M*(Up.*conj(Up));

if rand > pr(1)

k=2;

else

k=1;

end

W=(1/sqrt(pr(k)))*M(k,:)'.*Up;


Top

RSA Example


clear

%choose two large primes and calculate product

p=811

q=613

n=p*q

pause

%choose an integer d coprime to (p-1)(q-1)

d=3211

pmqm=(p-1)*(q-1)

gcd(d,pmqm)

pause

%Compute e the modular inverse

fp=factor(p-1)

fq=factor(q-1)

pause

phi=fix(pmqm*(1-1/2)*(1-1/3)*(1-1/5)*(1-1/17))

pause

e=modp(3211,phi-1,pmqm)

pause

%Broadcast e and n

% M is the integer message

M=123456

pause

%Encrypt the message

E=modp(M,e,n)

pause

%Send and then receiver decodes

Mdc=modp(E,d,n)


Top

s^p mod n  Function


function z=modp(s,p,a)

%Calculates s^p mod a for p < 2^20 ~ 10^6

c=zeros(1,20); c(1)=s;

for i=2:20

c(i)=mod(c(i-1)*c(i-1),a);

end

b=zeros(1,20);

for i=19:-1:0

r=p-2^i;

if r>=0

b(i+1)=1; p=r;

end

end

z=1;

for i=1:20

if b(i)==1

z=mod(z*c(i),a);

end

end


Top

Shor's Algorithm Simulation


% Shor's Algorithm on a conventional computer

%N = number to be factored

N=35

%choose a number relatively prime to N

x=23

gcd(x,N)

%choose a number ~ N^2

q=1024

a=1:q; fna=zeros(1,q);

fna(1)=x;

for i=2:(q)

fna(i)=mod(x*fna(i-1),N);

end

plot(a(1:25),fna(1:25))

grid on

pause

ffna=abs(ifft(fna));

plot(a,ffna)

axis([0,q,0,1])

ffna(1:500)

f=input('Enter Frequency ')

[Q,NU,D]=cfrac(f,1024);

Q,[NU;D]

[Q,NU,D]=cfrac(f+1,1024);

Q,[NU;D]

[Q,NU,D]=cfrac(f+2,1024);

Q,[NU;D]

[Q,NU,D]=cfrac(f+3,1024);

Q,[NU;D]

r=input('Enter Period ')

gcd(x^(r/2)-1,N)

gcd(x^(r/2)+1,N)


Top

Continued Fraction Expansion Function


function [Q,N,D] = cfrac(p,q)

%Calculates the Continued Fraction Expansion of p/q assuming that p<q

r=1; i=0;

while r>0

i=i+1;

Q(i)=floor(q/p);

r=q-Q(i)*p;

q=p; p=r;

end

n=length(Q);

for i=n:-1:2

a=Q(i); b=1;

for j=i-1:-1:1

c=a;

a=Q(j)*a+b;

b=c;

end

N(i)=b; D(i)=a;

end

N(1)=1; D(1)=Q(1);

end


Top