Louis Lamarche

Ellipticaly polarized electromagnetic wa...

Posted: Louis Lamarche 105 Product: Maple

September 14 2020

5 0

EllipticalyPolarizedWaveRefectionByMirrorFinal.mw

Typesetting vector names with up arrows ...

Posted: Louis Lamarche 105 Product: Maple 2018

January 10 2020

2 3

I was trying to display a Physics[Vectors] vector name in a 3dplot with an up arrow
on it. I found that this old 2008 trick still works in MAPLE 2018.

>

restart;

>	with(plots): with(Physics[Vectors]):

>	# Using MAPLE 2018.2

>	a:=arrow([-1,1,1],view=[-1.5..1.5,-1.5..1.5,-1.5..1.5]):

>	v_; t:= textplot3d([-1.1,1.1,1,v_]): display(a,t);

>	# I found this on an old 2008 post t:= textplot3d([-1.1,1.1,1,typeset(`#mover(mi(` \|\| v \|\| `),mo("→"))`)]): display(a,t);

Download VectorTypeSetting.mw

Quantum Teleportation in Dirac notation...

Posted: Louis Lamarche 105 Product: Maple

March 08 2018

4 4

This is an application of the previous posts
https://www.mapleprimes.com/posts/209057-Procedure-For-Expanding-Tensor-Product

I have a fourth version of the ExpandQop that will expand automaticaly the power of
quantum tensor product. This is just a minor change to the procedure.

Now here is an application for all this that will help understanding a little about
quantum computing. This is the classical concept of quantum teleportation.

You will need to run the above mentionned file and uncomment the save line in the file
before running the example.

LL

>

######################################################################
# NOTICE                                                             #
# Author: Louis Lamarche                                             #
#         Institute of Research of Hydro-Quebec (IREQ)               #
#         Science des données et haute performance                   #
#         2018, March 7                                              #
#                                                                    #
# Function name: ExpandQop (x)                                       #
#       Purpose: Compute the tensor product of two quantum           #
#                operators in Dirac notations                        #
#      Argument: x: a quantum operator                               #
# Improvements: Manage all +, -, *, /, ^, mod operations           #
#                in the argument. Manages multiple tensor products   #
#                like A*B*C*F                                        #
#       Version: 3                                                   #
#                                                                    #
# Copyrigth(c) Hydro-Quebec.                                        #
#        Note 1: Permission to use this softwate is granted if you   #
#                acknowledge its author and copyright                #
#        Note 2: Permission to copy this softwate is granted if you #
#                leave this 21 lines notice intact. Thank you.       #
######################################################################
restart;

>	with(Physics): interface(imaginaryunit=i): Setup(mathematicalnotation=true);

(1)

>	Setup(unitaryoperators={I,U,X,Y,Z,H,HI,CNOT,CnotI}); Setup(noncommutativeprefix={q,beta,psi});

(2)

>	Setup(bracketrules= { %Bracket(%Bra(q0), %Ket(q0))=1, %Bracket(%Bra(q1), %Ket(q1))=1, %Bracket(%Bra(q1), %Ket(q0))=0, %Bracket(%Bra(q0), %Ket(q1))=0 });

(3)

>

####################################################################################
# Load the procedure and set the required global variables
#
read "ExpandQop.m": optp:=op(0,Ket(q0)*Ket(q1)): optpx:= op(0,(Ket(q0)+Ket(q1))^2):
#
####################################################################################

>

#
# Pauli operators
#
print("Pauli gates");
I:=Ket(q0)*Bra(q0)+Ket(q1)*Bra(q1);        # = sigma[0]
X:=Ket(q1)*Bra(q0)+Ket(q0)*Bra(q1);        # = sigma[1] = sigma[x]
Y:=-i*Ket(q1)*Bra(q0)+i*Ket(q0)*Bra(q1);   # = sigma[2] = sigma[y]
Z:=Ket(q0)*Bra(q0)-Ket(q1)*Bra(q1);        # = sigma[3] = sigma[z]

(4)

>	############################## # Defining the Hadamard gate # ############################## print("Hadamard gate"); H:= Ket(q0)Bra(q0)/sqrt(2)+Ket(q0)Bra(q1)/sqrt(2)+Ket(q1)Bra(q0)/sqrt(2)-Ket(q1)Bra(q1)/sqrt(2);

(5)

>

# This is usefull to represent a 2 qubits system
# A more general approach is needed for a n qubit system.
DefineStates:=proc()
    Ket(q00):=Ket(q0)*Ket(q0); Ket(q01):=Ket(q0)*Ket(q1);
    Ket(q10):=Ket(q1)*Ket(q0); Ket(q11):=Ket(q1)*Ket(q1);
    Bra(q00):=Dagger(Ket(q00)); Bra(q01):=Dagger(Ket(q01));
    Bra(q10):=Dagger(Ket(q10)); Bra(q11):=Dagger(Ket(q11));
    return;
    end proc:
UndefineStates:=proc()
    Ket(q00):='Ket(q00)'; Ket(q01):='Ket(q01)';
    Ket(q10):='Ket(q10)'; Ket(q11):='Ket(q11)';
    Bra(q00):='Bra(q00)'; Bra(q01):='Bra(q01)';
    Bra(q10):='Bra(q10)'; Bra(q11):='Bra(q11)';
    return;
    end proc:

>	#################################### # Defining the CNOT gate (2 qubits) #################################### print("CNOT gate"); CNOT:=Ket(q00)Bra(q00)+ Ket(q01)Bra(q01)+ Ket(q11)Bra(q10)+Ket(q10)Bra(q11); DefineStates(); 'CNOT'=CNOT;

(6)

>

###########################
# Defining the Bell states
###########################
Ket(beta,x,y)='CNOT.(((H.Ket(x)))*Ket(y))';
Ket(beta00):=CNOT.(Expand((H.Ket(q0)))*Ket(q0));
Ket(beta01):=CNOT.(Expand((H.Ket(q0)))*Ket(q1));
Ket(beta10):=CNOT.(Expand((H.Ket(q1)))*Ket(q0));
Ket(beta11):=CNOT.(Expand((H.Ket(q1)))*Ket(q1));

(7)

>

##########################################################
# Quantum teleportation
# Reference: Quantum Computation and Quantum Information
#            10th Anniversary Edition
#            Michael A. Nielsen & Isaac L. Chuang
#            Cambridge University Press, Cambridge 2010
#            pp 25-28
##########################################################
print("State to be teleported");
Ket(psi) := a*Ket(q0)+b*Ket(q1);
print("Step 1: Compute the tensor product of the state to be teleported with ", 'Ket(beta00)');
Ket(psi[0])='Ket(psi)'*'Ket(beta00)';
Ket(psi[0]):=Expand(Ket(psi)*Ket(beta00));
print("This is a 3 qubits state");
#######
print("Step 2: Pass these 3 qubits through a CNOT*I operator");
'CnotI'='CNOT*I';
CnotI:=ExpandQop(Expand(CNOT*I)):
#
# To see what the CNOTI operator looks like
#
# print("CNOTI=");
# print(op(1,CNOTI)+op(2,CNOTI)+op(3,CNOTI)+op(4,CNOTI));
# print(op(5,CNOTI)+op(6,CNOTI)+op(7,CNOTI)+op(8,CNOTI));
'Ket(psi[1])'='CnotI.Ket(psi[0])';
Ket(psi[1]):=Expand(CnotI.Ket(psi[0]));
#######
print("Step 3: Pass these 3 qubits through an Haldamard*I operator");
'HalI'='H*I';
HalI:=ExpandQop(Expand(H*I)):
#
# To see what the Haldamard*I operator looks like
#
# print("HalI=");
# print(op(1,HalI)+op(2,HalI)+op(3,HalI)+op(4,HalI));
# print(op(5,HalI)+op(6,HalI)+op(7,HalI)+op(8,HalI));
'Ket(psi[2])'='HalI.Ket(psi[1])';
Ket(psi[2]):=Expand(HalI.Ket(psi[1]));

(8)

>

UndefineStates();
print("Using contracted names for the first two qubits");
Ket(q00)*Bra(q0)*Bra(q0)='I';
Ket(q01)*Bra(q0)*Bra(q1)='I';
Ket(q10)*Bra(q1)*Bra(q0)='I';
Ket(q11)*Bra(q1)*Bra(q1)='I';
'Ket(psi[2])'=Ket(q00)*Bra(q0)*Bra(q0).Ket(psi[2])+
              Ket(q01)*Bra(q0)*Bra(q1).Ket(psi[2])+
              Ket(q10)*Bra(q1)*Bra(q0).Ket(psi[2])+
              Ket(q11)*Bra(q1)*Bra(q1).Ket(psi[2]);

(9)

>	print("Rewriting this result by hand"); 'Ket(psi[2])'=(Ket(q00)(aKet(q0)+bKet(q1))+ Ket(q01)(aKet(q0)-bKet(q1))+ Ket(q10)(aKet(q1)+bKet(q0))+ Ket(q11)(aKet(q1)-bKet(q0)))/2;

Physics:-Ket(psi[2]) = (1/2)*Physics:-`*`(Physics:-Ket(q00), a*Physics:-Ket(q0)+b*Physics:-Ket(q1))+(1/2)*Physics:-`*`(Physics:-Ket(q01), a*Physics:-Ket(q0)-b*Physics:-Ket(q1))+(1/2)*Physics:-`*`(Physics:-Ket(q10), a*Physics:-Ket(q1)+b*Physics:-Ket(q0))+(1/2)*Physics:-`*`(Physics:-Ket(q11), a*Physics:-Ket(q1)-b*Physics:-Ket(q0))

(10)

>

DefineStates();
print("If Alice measures 00 Bob does noting");
''I'.   '2*Bra(q00).Ket(psi[2])'' = I.   2*Bra(q00).Ket(psi[2]);
print("If Alice measures 01 Bob applies the X gate");
''X'.   '2*Bra(q01).Ket(psi[2])'' = X.   2*Bra(q01).Ket(psi[2]);
print("If Alice measures 10 Bob applies the Z gate");
''Z'.   '2*Bra(q10).Ket(psi[2])'' = Z.   2*Bra(q10).Ket(psi[2]);
print("If Alice measures 11 Bob applies the X gate and then the Z gate");
''Z'.'X'. '2*Bra(q11).Ket(psi[2])'' = Z.X. 2*Bra(q11).Ket(psi[2]);

(11)

>

Download QuantumTeleportation.mw

Procedure for expanding tensor product o...

Posted: Louis Lamarche 105 Product: Maple

March 08 2018

1 1

Version 2 do not enable to expand multiple product like A*A*B*E
Version 3 will now do that.
I just forgot to add this feature.

LL.

>

######################################################################
# NOTICE                                                             #
# Author: Louis Lamarche                                             #
#         Institute of Research of Hydro-Quebec (IREQ)               #
#         Science des données et haute performance                   #
#         2018, March 7                                              #
#                                                                    #
# Function name: ExpandQop (x)                                       #
#       Purpose: Compute the tensor product of two quantum           #
#                operators in Dirac notations                        #
#      Argument: x: a quantum operator                               #
# Improvements: Manage all +, -, *, /, ^, mod operations           #
#                in the argument. Manages multiple tensor products   #
#                like A*B*C*F                                        #
#       Version: 3                                                   #
#                                                                    #
# Copyrigth(c) Hydro-Quebec.                                        #
#        Note 1: Permission to use this softwate is granted if you   #
#                acknowledge its author and copyright                #
#        Note 2: Permission to copy this softwate is granted if you #
#                leave this 21 lines notice intact. Thank you.       #
######################################################################
restart;

>	with(Physics): interface(imaginaryunit=i): Setup(mathematicalnotation=true);

(1)

>	Setup(quantumoperators={A,B,C,Cn}); Setup(noncommutativeprefix={a,b,q});

(2)

>

opexp:= op(0,10^x):            # exponentiation id
opnp := op(0,10*x):            # normal product id
optp := op(0,Ket(q0)*Ket(q1)): # tensor product id
opdiv:= `Fraction`:            # fraction       id
opsym:= op(0,x):               # symbol         id
opint:= op(0,10):              # integer        id
opflt:= op(0,10.0):            # float          id
opcpx:= op(0,i):               # complex        id
opbra:= op(0,Bra(q)):          # bra            id
opket:= op(0,Ket(q)):          # ket            id
opmod:= op(0, a mod b):        # mod            id
ExpandQop:=proc(x)
    local nx,ret,j,lkb,cbk,rkb,no,lop,success;
    lop:=op(0,x);
    no:=nops(x);
    if lop = opsym or lop = opint or lop = opflt or
       lop = opbra or lop = opket or lop = opcpx then
         ret:=x;
    else
    if lop = opexp then
        ret:=x;
    else
    if lop = opnp then
        ret:=1;
        for j from 1 to no do
            ret:=ret*ExpandQop(op(j,x));
        end do;
    else
    if lop = `+` then
        ret:=0;
        for j from 1 to no do
            ret:=ret+ExpandQop(op(j,x));
        end do;
    else
    if lop = `-` then
        ret:=0;
        for j from 1 to no do
            ret:=ret-ExpandQop(op(j,x));
        end do;
    else
    if lop = opdiv then
       ret:=1;
       for j from 1 to no do
           ret:=ret/ExpandQop(op(j,x));
       end do;
    else
    if lop = opmod then
       ret:=x;
    else
    if lop = optp then
       if (no > 3 ) then
           success:=false;
           nx:=x;
           while not success do
             lkb:=0; cbk:=0; rkb:=0;ret:=1;
             for j from 1 to no do
                 if (j>1) then
                      if(lkb=0) then
                          if( type(op(j-1,nx),Ket) and
                              type(op(j,nx),Bra) ) then lkb:=j-1; fi;
                      else
                          if( type(op(j-1,nx),Ket) and
                              type(op(j,nx),Bra) ) then rkb:=j;   fi;
                      fi;
                      if( type(op(j-1,nx),Bra) and type(op(j,nx),Ket) )
                                                   then cbk:=j;   fi;
                 fi;
             end do;
             if ( (lkb < cbk) and (cbk<rkb) ) then
                 for j from 1     to lkb   do ret := ret*op(j,nx); end do;
                 for j from cbk   to no    do ret := ret*op(j,nx); end do;
                 for j from lkb+1 to cbk-1 do ret := ret*op(j,nx); end do;
             else
               ret:=nx;
             fi;

             if nx = ret then
                success := true;
             else
                nx := ret;
             fi
           end do;
       else
           ret:=x;
       fi;
    else ret:=x;
    fi; # optp
    fi; # opmod
    fi; # opdiv
    fi; # `-`
    fi; # `+`
    fi; # `opnp
    fi; # `opexp`
    fi; # opsym, opint, opflt, opbra, opket, opcpx

    return ret;
end proc:

# For saving
# save opexp,opnp,optp,opdiv,opint,opflt,opcpx,opbra,opket,opmod, ExpandQop,"ExpandQop.m"

>

# Let A be an operator in a first Hilbert space of dimension n
# using the associated orthonormal ket and bra vectors
#
#
kets1:=Ket(a1)*Ket(a2)*Ket(a3)*Ket(a4)*Ket(a5):
A:=kets1*Dagger(kets1);

# Let B be an operator in a second Hilbert (Ketspace of dimension m
# using the associated orthonormal ket and bra vectors
#
#
kets2:=Ket(b1)*Ket(b2)*Ket(b3):
B:=kets2*Dagger(kets2);

# The tensor product of the two operators acts on a n+m third
# Hilbert space unsing the appropriately ordered ket
# and bra vectors of the two preceding spaces. The rule for
# building this operator in Dirac notation is as follows,
#
#

print("Maple do not compute the tensor product of operators,");
print("C=A*B gives:");
C:=A*B;

print("ExpandQop(C) gives the expected result:");
Cn:=ExpandQop(C);

(3)

>	kets3:=kets1*kets2; bras3:=Dagger(kets3); print("Matrix elements computed with C appears curious"); 'bras3.C. kets3'="..."; bras3.C.kets3; print("Matrix elements computed with Cn as expected"); 'bras3.Cn.kets3'=bras3.Cn.kets3;

(4)

>	print("Example"); En:=ExpandQop(10(1-x+y+z)i(1/sqrt(2))A*B);

(5)

>	print("Another example"); 'F'='AB/sqrt(2)+BA/sqrt(2)'; F:=AB/sqrt(2)+BA/sqrt(2): 'op(1,F)'=op(1,F); 'op(2,F)'=op(2,F); 'Fn'='ExpandQop(F)'; Fn:=ExpandQop(F): 'op(1,Fn)'=op(1,Fn); 'op(2,Fn)'=op(2,Fn);

(6)

>	print("Final example, multiple products"); G:=BBB; 'G'=ExpandQop(G);

(7)

>

Download ExpandQopV3.mw

Procedure for expanding tensor product o...

Posted: Louis Lamarche 105 Product: Maple

March 07 2018

2 8

Here is a major upgrade of the procedure I submitted in february.

https://www.mapleprimes.com/posts/209030-Procedure-For-Computing-The-Tensor-Product

There is a line after the procedure to save it in the file "ExpandQop.m"
In future post I will use it in order to minimize the size of the examples.

Louis Lamarche

>

######################################################################
# NOTICE                                                             #
# Author: Louis Lamarche                                             #
#         Institute of Research of Hydro-Quebec (IREQ)               #
#         Science des données et haute performance                   #
#         2018, March 7                                              #
#                                                                    #
# Function name: ExpandQop (x)                                       #
#       Purpose: Compute the tensor product of two quantum           #
#                operators in Dirac notations                        #
#      Argument: x: a quantum operator                               #
# Improvements: Manage all +, -, *, /, ^, mod operations           #
#                in the argument                                     #
#       Version: 2                                                   #
#                                                                    #
# Copyrigth(c) Hydro-Quebec.                                        #
#        Note 1: Permission to use this softwate is granted if you   #
#                acknowledge its author and copyright                #
#        Note 2: Permission to copy this softwate is granted if you #
#                leave this 21 lines notice intact. Thank you.       #
######################################################################
restart;

>	with(Physics): interface(imaginaryunit=i): Setup(mathematicalnotation=true);

(1)

>	Setup(quantumoperators={A,B,C,Cn}); Setup(noncommutativeprefix={a,b,q});

(2)

>

opexp:= op(0,10^x):            # exponentiation id
opnp := op(0,10*x):            # normal product id
optp := op(0,Ket(q0)*Ket(q1)): # tensor product id
opdiv:= `Fraction`:            # fraction       id
opsym:= op(0,x):               # symbol         id
opint:= op(0,10):              # integer        id
opflt:= op(0,10.0):            # float          id
opcpx:= op(0,i):               # complex        id
opbra:= op(0,Bra(q)):          # bra            id
opket:= op(0,Ket(q)):          # ket            id
opmod:= op(0, a mod b):        # mod            id
ExpandQop:=proc(x)
    local ret,j,lkb,cbk,rkb,no,lop;
    lkb:=0; cbk:=0; rkb:=0;
    lop:=op(0,x);
    no:=nops(x);
    if lop = opsym or lop = opint or lop = opflt or
       lop = opbra or lop = opket or lop = opcpx then
         ret:=x;
    else
    if lop = opexp then
        ret:=x;
    else
    if lop = opnp then
        ret:=1;
        for j from 1 to no do
            ret:=ret*ExpandQop(op(j,x));
        end do;
    else
    if lop = `+` then
        ret:=0;
        for j from 1 to no do
            ret:=ret+ExpandQop(op(j,x));
        end do;
    else
    if lop = `-` then
        ret:=0;
        for j from 1 to no do
            ret:=ret-ExpandQop(op(j,x));
        end do;
    else
    if lop = opdiv then
       ret:=1;
       for j from 1 to no do
           ret:=ret/ExpandQop(op(j,x));
       end do;
    else
    if lop = opmod then
       ret:=x;
    else
    if lop = optp then
        ret:=1;
       if (no > 3 ) then
           for j from 1 to no do
               if (j>1) then
                    if(lkb=0) then
                        if( type(op(j-1,x),Ket) and
                            type(op(j,x),Bra) ) then lkb:=j-1; fi;
                    else
                        if( type(op(j-1,x),Ket) and
                            type(op(j,x),Bra) ) then rkb:=j;   fi;
                    fi;
                    if( type(op(j-1,x),Bra) and type(op(j,x),Ket) )
                                                then cbk:=j;   fi;
               fi;
           end do;
           if ( (lkb < cbk) and (cbk<rkb) ) then
               for j from 1     to lkb   do ret := ret*op(j,x); end do;
               for j from cbk   to no    do ret := ret*op(j,x); end do;
               for j from lkb+1 to cbk-1 do ret := ret*op(j,x); end do;
           else
               ret:=x;
           fi;
       else
           ret:=x;
       fi;
    else ret:=x;
    fi; # optp
    fi; # opmod
    fi; # opdiv
    fi; # `-`
    fi; # `+`
    fi; # `opnp
    fi; # `opexp`
    fi; # opsym, opint, opflt, opbra, opket, opcpx

    return ret;
end proc:

# For saving
# save opexp,opnp,optp,opdiv,opint,opflt,opcpx,opbra,opket,opmod, ExpandQop,"ExpandQop.m"

>

# Let A be an operator in a first Hilbert space of dimension n
# using the associated orthonormal ket and bra vectors
#
#
kets1:=Ket(a1)*Ket(a2)*Ket(a3)*Ket(a4)*Ket(a5):
A:=kets1*Dagger(kets1);

# Let B be an operator in a second Hilbert (Ketspace of dimension m
# using the associated orthonormal ket and bra vectors
#
#
kets2:=Ket(b1)*Ket(b2)*Ket(b3):
B:=kets2*Dagger(kets2);

# The tensor product of the two operators acts on a n+m third
# Hilbert space unsing the appropriately ordered ket
# and bra vectors of the two preceding spaces. The rule for
# building this operator in Dirac notation is as follows,
#
#

print("Maple do not compute the tensor product of operators,");
print("C=A*B gives:");
C:=A*B;

print("ExpandQop(C) gives the expected result:");
Cn:=ExpandQop(C);

(3)

>	kets3:=kets1*kets2; bras3:=Dagger(kets3); print("Matrix elements computed with C appears curious"); 'bras3.C. kets3'="..."; bras3.C.kets3; print("Matrix elements computed with Cn as expected"); 'bras3.Cn.kets3'=bras3.Cn.kets3;

(4)

>	print("Example"); En:=ExpandQop(10(1-x+y+z)i(1/sqrt(2))A*B);

(5)

>	print("Another example"); 'F'='AB/sqrt(2)+BA/sqrt(2)'; F:=AB/sqrt(2)+BA/sqrt(2): 'op(1,F)'=op(1,F); 'op(2,F)'=op(2,F); 'Fn'='ExpandQop(F)'; Fn:=ExpandQop(F): 'op(1,Fn)'=op(1,Fn); 'op(2,Fn)'=op(2,Fn);