# OK to post homework
# James Betti, 23 Feb 2025, Assignment 8

for _ to 3 do
    M := RHM(3,10);
    V := RSV(3,10);
    print(evalf(ProbDist(M,V)[3]), (evalf(ExpMV(M,V)))); # should be close
od:

for _ to 3 do
    V := RSV(2,100);
    print(add(ExpMV(PM()[i],V)^2,i=1..3)); # should be 1
od:

EntanglementAlice := proc(V) local A,i;
    A := map(P->TP(P,PM()[4]),PM());
    add(ExpMV(A[i],V)^2,i=1..3);
end:

RandomProductState4 := proc(K);
    op(TP([RSV(2,K)],[RSV(2,K)]));
end:

[seq(EntanglementAlice(ES()[i]),i=1..4)];
evalf(EntanglementAlice(RandomProductState4(10)));

maxV := [0,[]]:
minV := [infinity,[]]:
for _ to 20 do
    V := RSV(4,10);
    x := EntanglementAlice(V);
    if x > maxV[1] then
        maxV := [x,V];
    fi;
    if x < minV[1] then
        minV := [x,V];
    fi;
od:
maxV, minV;                             # 0.9573, 0.0873