#Feb. 10, 2025 C6.txt Help6:=proc(): print(`Mul(A,B), Com(A,B), RSV(n,K), Del(L,A), CheckUP(A,B,V) `): end: #Mul(A,B): the product of matrix A and B (assuming that it exists) Mul:=proc(A,B) local i,j,k,n: n:=nops(A): [seq([seq( add(A[i][k]*B[k][j] , k=1..n) , j=1..n)],i=1..n)]: end: #Com(A,B: The commutator of A and B Com:=proc(A,B): Mul(A,B)-Mul(B,A):end: #RSV(n,K): RSV:=proc(n,K) local ra,i,v,a: ra:=rand(-K..K): v:=[seq(ra()+I*ra(),i=1..n)]: a:=sqrt(add(v[i]*conjugate(v[i]),i=1..n)): [seq(v[i]/a,i=1..n)]: end: #Del(L,A): the "avergae error", i.e. standard deviation pf the prob. distribution of #measuring the observable L when the state is A Del:=proc(L,A) local n,pd,ev,i,av: n:=nops(L): pd:=ProbDist(L,A): ev:=pd[1]: av:=pd[3]: pd:=pd[2]: sqrt(add(ev[i]^2*pd[i],i=1..n)-av^2): end: #CheckUP(A,B,V): the ratio of Del(A,V)*Del(B,V)/Expectation(Com(A,B),V) (it should be >=1/2) CheckUP:=proc(A,B,V) local C: C:=Com(A,B): if abs(ProbDist(C,V)[3])=0 then RETURN(FAIL): fi: evalf(Del(A,V)*Del(B,V)/abs(ProbDist(C,V)[3])): end: #old stuff #C5.txt, Feb. 6, 2025 Experimental Mathematics (Rutgers, Dr. Z.) Help5:=proc(): print(` PM(), RHM(n,K) , EV(L), EVC1(L,g), EVC(L), IP(A,B), ProbDist(L,A) `): end: with(linalg): #PM(): the three famous Pauli matrices sigma_x,sigma_y, sigma_z PM:=proc(): [ [[0,1],[1,0]], [[0,-I],[I,0]], [[1,0],[0,-1]] ]: end: #RHM(n,K): a random n by n Hermitian matrix with entriries from -K to K RHM:=proc(n,K) local ra,i,j,L: ra:=rand(-K..K): for i from 1 to n do for j from 1 to i-1 do L[i,j]:=ra()+I*ra(): L[j,i]:=conjugate(L[i,j]): od: L[i,i]:=ra(): od: [seq([seq(L[i,j],j=1..n)],i=1..n)]: end: #EV(L): inputs a Hermitian matrix (say n by n) and outputs the LIST of eigenvalues EV:=proc(L) local x,L1,i,n: n:=nops(L): L1:=[seq(L[i]-[0$(i-1),x,0$(n-i)],i=1..n)]: [solve(det(L1),x)]: end: #EVC1(L,g): inputs a matrix L and an eigenvalue g finds the NORMALIZED eigenvector EVC1:=proc(L,g) local v,i,j,var,eq,V,n,L1,V1,va,k: n:=nops(L): L1:=[seq(L[i]-[0$(i-1),g,0$(n-i)],i=1..n)]: if simplify(det(L1))<>0 then RETURN(FAIL): fi: V:=[seq(v[i],i=1..n)]: var:={seq(v[i],i=1..n)}: eq:={seq(add(L1[i,j]*V[j],j=1..n),i=1..n)}: var:=solve(eq,var): #V1 the set of unknowns that is arbitrary V1:={}: for va in var do if op(1,va)=op(2,va) then V1:=V1 union {op(1,va)}: fi: od: if nops(V1)<>1 then RETURN(FAIL): fi: V:=subs(var,V): V:=subs(V1[1]=1,V): k:=sqrt(add(abs(V[i])^2,i=1..n)): [seq(V[i]/k,i=1..n)]: end: #EVC(L): inputs a Hermitian matrix and outputs the pair ListOfEigenvalues, ListOfEigenvectors EVC:=proc(L) local n,ev,evc,i: n:=nops(L): ev:=EV(L): evc:=[seq(EVC1(L,ev[i]),i=1..n)]: [ev,evc]: end: #IP(A,B): the inner product of state-vectors A and B IP:=proc(A,B) local n,i: n:=nops(A): if nops(B)<>n then RETURN(FAIL): fi: simplify(add(A[i]*conjugate(B[i]),i=1..n)): end: #ProbDist(L,A): Given an obervalble represented by the Hermitian matrix L and a state vector A #outputs first all the eigenvectors (pure states) and the prob. of winding in the respective pure states #followed by the expected value of A (i.e. if you make the observation L and the state is A the average in the long run ProbDist:=proc(L,A) local evc,i,A1,k,PD,n: n:=nops(L): k:=sqrt(add(abs(A[i])^2,i=1..n)): A1:=[seq(A[i]/k,i=1..n)]: evc:=EVC(L): PD:=[seq(abs(IP(A1,evc[2][i]))^2,i=1..n)]: [evc[1],PD,add(PD[i]*evc[1][i],i=1..n)]: end: