Hello,
I just started to use qetlab but whtever I tried I always an error. For intance:
When I want to try that:

GHZState(4,7,[1,2,3,4]/sqrt(30))

I have this error:

opt_args' is not found in the current folder or on the MATLAB path, but exists in:
/MATLAB Drive/QETLAB-0.9/helpers

Change the MATLAB current folder or add its folder to the MATLAB path.

Error in GHZState (line 36)
[coeff] = opt_args({ ones(1,dim)/sqrt(dim) },varargin{:});

However When I try the same comment under the helper this time I have that:

GHZState' is not found in the current folder or on the MATLAB path, but exists in:
/MATLAB Drive/QETLAB-0.9

Change the MATLAB current folder or add its folder to the MATLAB path.

What I mean is that wherever I try I have an error

Or let's take a simple example:
rho = RandomDensityMatrix(9);
This time I have that:

'opt_args' is not found in the current folder or on the MATLAB path, but exists in:
/MATLAB Drive/QETLAB-0.9/helpers

Change the MATLAB current folder or add its folder to the MATLAB path.

Error in RandomDensityMatrix (line 33)
[re,k,dist] = opt_args({ 0, dim, 'haar' },varargin{:});

First of all, thanks a lot for documenting LocalDistinguishability on the wiki.

In an ideal workflow, documentation should go hand by hand with the code, that is, there should be documentation changes in each Pull Request where you change the signature of the functions in the code. I am opening this issue to discuss a way to automatize documentation.

As we discussed by email, this is not a priority, but I believe it's still good to write down potential solutions:

• Mediawiki as it is now (source on Github)
• MATLAB documentation hosted on Github Pages
  • free hosting
• readme.io
  • free for open source projects

(more pro and cons to be added)

In PartialTranspose.m there is a check on the input matrix to see if it is numeric or not:

if(isnumeric(Xpt)) % if the input is a numeric matrix, perform the partial transpose operation the fastest way we know how
    Xpt = cell2mat(mat2cell(Xpt, sub_sys_vecR, sub_sys_vecC).'); % partial transpose on first subsystem
else % if the input is not numeric (such as a variable in a semidefinite program), do a slower method that avoids mat2cell (mat2cell doesn't like non-numeric arrays)
    Xpt = reshape(permute(reshape(Xpt,[sub_sys_vecR(1),sub_prodR,sub_sys_vecC(1),sub_prodC]),[1,4,3,2]),[prod_dimR,prod_dimC]);
end

Consider the following (numeric) matrix:

X = reshape(1:210,14,15)'

When we encounter the above if statement, the first condition is triggered and we obtain the following partial transposed matrix

     1     2     3     4     5     6     7    71    72    73    74    75    76    77   141   142   143   144   145   146   147
    15    16    17    18    19    20    21    85    86    87    88    89    90    91   155   156   157   158   159   160   161
    29    30    31    32    33    34    35    99   100   101   102   103   104   105   169   170   171   172   173   174   175
    43    44    45    46    47    48    49   113   114   115   116   117   118   119   183   184   185   186   187   188   189
    57    58    59    60    61    62    63   127   128   129   130   131   132   133   197   198   199   200   201   202   203
     8     9    10    11    12    13    14    78    79    80    81    82    83    84   148   149   150   151   152   153   154
    22    23    24    25    26    27    28    92    93    94    95    96    97    98   162   163   164   165   166   167   168
    36    37    38    39    40    41    42   106   107   108   109   110   111   112   176   177   178   179   180   181   182
    50    51    52    53    54    55    56   120   121   122   123   124   125   126   190   191   192   193   194   195   196
    64    65    66    67    68    69    70   134   135   136   137   138   139   140   204   205   206   207   208   209   210

However, if we comment out the if condition and instead trigger the else condition, then the resulting matrix is

     1     9     4    12     7    78    73    81    76    84   142   150   145   153
    15    23    18    26    21    92    87    95    90    98   156   164   159   167
    29    37    32    40    35   106   101   109   104   112   170   178   173   181
    43    51    46    54    49   120   115   123   118   126   184   192   187   195
    57    65    60    68    63   134   129   137   132   140   198   206   201   209
     8     3    11     6    14    72    80    75    83   141   149   144   152   147
    22    17    25    20    28    86    94    89    97   155   163   158   166   161
    36    31    39    34    42   100   108   103   111   169   177   172   180   175
    50    45    53    48    56   114   122   117   125   183   191   186   194   189
    64    59    67    62    70   128   136   131   139   197   205   200   208   203
     2    10     5    13    71    79    74    82    77   148   143   151   146   154
    16    24    19    27    85    93    88    96    91   162   157   165   160   168
    30    38    33    41    99   107   102   110   105   176   171   179   174   182
    44    52    47    55   113   121   116   124   119   190   185   193   188   196
    58    66    61    69   127   135   130   138   133   204   199   207   202   210

One would assume that these two conditions should yield the same matrix of the same dimension. Perhaps I have a misunderstanding about the way in which this is working?

Feel free to close this issue if my assumption here is wrong, but I believe in Tensor.m, if you create a vector say in the following way:

q0 = [1;0]; Tensor(q0,0) Tensor(q0,1)

Then I would imagine that the first Tensor call would return a scalar or "1", while the second would return the vector "q0". The second call to the Tensor function works as expected, but the first also returns the vector "q0" as opposed to just "1".

Again, if my intuition is incorrect here, please feel free to close this issue. Thanks!

It seems as if the ChoiMatrix function may be performing a slightly off conversion for certain cases when the input to the function is a set of Kraus operators.

For instance:

% The Choi matrix of the depolarizing channel on two-dimensions.
full(DepolarizingChannel(2))
ans =

    0.5000         0         0         0
         0    0.5000         0         0
         0         0    0.5000         0
         0         0         0    0.5000

% Obtain the Kraus operators corresponding to the depolarizing channel.
kraus_ops = KrausOperators(DepolarizingChannel(2));

T = {kraus_ops{1}, kraus_ops{1}';  
       kraus_ops{2}, kraus_ops{2}';
       kraus_ops{3}, kraus_ops{3}';
       kraus_ops{4}, kraus_ops{4}';
};

% Compute the Choi matrix corresponding to the Kraus operators. I *believe* we should
% expect to get the Choi matrix of the depolarizing channel back, but instead, we get
% the following matrix.
ChoiMatrix(T)
ans =

    0.5000         0         0         0
         0         0    0.5000         0
         0    0.5000         0         0
         0         0         0    0.5000

I might be misinterpreting/misusing something, but just in case I am not, I figured I would mention that here.

I'll start this issue off by saying there is a high likelihood that there is some form of user error on my part happening here.

To give context, I am attempting to verify that the operators from Equation (35) in arXiv:1408.6981 are block-positive via QETLAB.

The following is a script that attempts to verify the above fact. The equation reference numbers present in the arXiv:1408.6981 posting

% Define variable epsilon \in [0,1] for resource state.
eps = 0.5;

% |0>, |1>.
e0 = [1;0]; e1 = [0;1];

% Define the standard two-qubit Bell states (Equation (2)).
psi_1 = [1/sqrt(2); 0; 0; 1/sqrt(2)]; 
psi_2 = [1/sqrt(2); 0; 0; -1/sqrt(2)];
psi_3 = [0; 1/sqrt(2); 1/sqrt(2); 0]; 
psi_4 = [0; 1/sqrt(2); -1/sqrt(2); 0];

% Define the corresponding density matrices of the Bell states.
rho_1 = psi_1 * psi_1'; rho_2 = psi_2 * psi_2';
rho_3 = psi_3 * psi_3'; rho_4 = psi_4 * psi_4';

% Define the resource state (Equation (18)):
eps_p = sqrt((1+eps)/2); eps_m = sqrt((1-eps)/2);
tau_vec = [eps_p; 0; 0; eps_m];
tau = tau_vec * tau_vec';

% Hermitian operator (Equation (32)).
H = 1/3 * (kron(eye(4), tau)/2 + ... 
    sqrt(1 - eps^2) * kron(rho_4, PartialTranspose(rho_4, 1, [2, 2])));

% Verify that the trace of H is equal to the closed-form expression (Equation (33)).
trace(H)
eq = 1/3 * (2 + sqrt(1-eps^2))

% Define Hermitian operators Q_{k,\eps} as in Equation (35) for k=1. 
% Should swap the second and third subsystems. 
% Each of the four subsystems is 2-dimensional.
Q0 = Swap(H - 1/3 * kron(rho_1, tau), [2, 3], [2, 2, 2, 2]);

% For k = 1,2,3, these operators should be block-positive however this 
% returns `0` which indicates that it is not determined to be block-positive.
IsBlockPositive(Q0)

I am assuming that perhaps the way in which I'm applying the Swap operator is perhaps not correct, although from my understanding I believe this is the correct way to go about doing this.

This might be related as well, but for eps=0.0 this appears to yield an Inaccurate state for the CVX solver. Could be related to my overall potential misunderstanding here as well though.

It's always a good idea to keep a change log.

Embarassingly enough the bug is in my own function :)
Check this script: https://gist.github.com/cosenal/65bfaec556e7d20855ea
It should output 5/6, instead it outputs 1 (states are from http://arxiv.org/abs/quant-ph/0103098)
It probably fails in many other cases too.

@nathanieljohnston @vprusso Can you help me with debugging this?

Looks like the following line in PermuteSystems.m:

Is checking for a different property than what is being stated in the comment below--that is, whether the length of perm and dim are equal. I think the correct fix here would be to alter the above line to:

if length(dim) ~= num_sys

Obviously, this won't cause any correctness issue, but it may incorrectly flag an invalid input as valid or vice-versa.

Just realized that the operator-Schmidt decomposition function OperatorSchmidtDecomposition works on the basis of arbitrary operators, but Hermitian ones.

This runs contrary to the definition e.g. in http://www.njohnston.ca/2014/06/what-the-operator-schmidt-decomposition-tells-us-about-entanglement/

One way to perform this decomposition is to first decompose the Hermitian operator on a basis of (local) Hermitian matrices such that

rho = \sum_{ij} c_{ij} A_i (x) B_j

where the c_{ij} are real, perform the SVD on c_{ij}, and then expand back the orthonormal Hermitian operators from the result [U S V] = svd(c). This is how I'm doing it for myself, but I'd like to submit a PR for QETLAB.

Is there a better way?

@nathanieljohnston It would be awesome if you could create a DOI for this repository. It should be very easy, just follow https://guides.github.com/activities/citable-code/

This way I could cite it properly in my thesis. I am pretty sure this needs to be done by the owner of the repository, otherwise I would have done it myself. Thanks!