Repository: qutip/qutip (python)

Issue:
Error occured in tensor_contract
### Bug Description

When attempting to discard the last system in a superoperator, an error occurs if the system sizes are unequal. I think this issue may be caused by an incorrect contraction of dimensions.

### Code to Reproduce the Bug

```shell
from qutip import *
A = qeye([[1], [1, 3]])
tensor_contract(A, (2, 5))
```


### Code Output

```shell
Traceback (most recent call last):
  File "from qutip import *.py", line 4, in <module>
    tensor_contract(A, (2, 5))
  File "/opt/anaconda3/envs/qutip/lib/python3.12/site-packages/qutip/core/tensor.py", line 341, in tensor_contract
    qmtx = qtens.reshape((l_mtx_dims, r_mtx_dims))
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
ValueError: cannot reshape array of size 9 into shape (1,1)
```


### Expected Behaviour

It should output a matrix with a shape (1, 1) instead of size 9. 

### Your Environment

```shell
QuTiP Version:      5.0.4
Numpy Version:      2.1.3
Scipy Version:      1.14.1
Cython Version:     None
Matplotlib Version: 3.9.2
Python Version:     3.12.7
Number of CPUs:     8
BLAS Info:          Generic
INTEL MKL Ext:      False
Platform Info:      Darwin (arm64)
```


### Additional Context

_No response_

Pull request description:
Fix bug in permutation of tensor reprensentation
**Description**
Fix a bug in the function to get the permutation to get a tensor from a uneven super operators.  (`[ [[2, 3], [4, 5, 6]], [...]]`) 

Add tests for `to_tensor_rep` and `from_tensor_rep`.

**Related issues or PRs**
fix #2555

Patch:
diff --git a/doc/changes/2555.bugfix b/doc/changes/2555.bugfix
new file mode 100644
index 000000000..dbcdfa36b
--- /dev/null
+++ b/doc/changes/2555.bugfix
@@ -0,0 +1,1 @@
+Fix getting tensor permutation for uneven super operators.
diff --git a/qutip/core/dimensions.py b/qutip/core/dimensions.py
index 56fe35396..3ef5bd37a 100644
--- a/qutip/core/dimensions.py
+++ b/qutip/core/dimensions.py
@@ -250,7 +250,7 @@ def dims_to_tensor_shape(dims):
     tensor_shape : tuple
         NumPy shape of the corresponding tensor.
     """
-    perm = dims_to_tensor_perm(dims)
+    perm = np.argsort(dims_to_tensor_perm(dims))
     dims = flatten(dims)
     return tuple(map(partial(getitem, dims), perm))
 
@@ -936,10 +936,10 @@ class Dimensions(metaclass=MetaDims):
         # dims_to_tensor_perm
         stepl = self.to_.step()
         stepr = self.from_.step()
-        return list(np.concatenate([
+        return list(np.argsort(np.concatenate([
             np.argsort(stepl)[::-1],
             np.argsort(stepr)[::-1] + len(stepl)
-        ]))
+        ])))
 
     def remove(self, idx: int | list[int]) -> "Dimensions":
         """


Test patch:
diff --git a/qutip/tests/core/test_dimensions.py b/qutip/tests/core/test_dimensions.py
index 96b2ed4cd..72d7aa246 100644
--- a/qutip/tests/core/test_dimensions.py
+++ b/qutip/tests/core/test_dimensions.py
@@ -6,7 +6,8 @@ import qutip
 from qutip.core.dimensions import (
     flatten, unflatten, enumerate_flat, deep_remove, deep_map,
     dims_idxs_to_tensor_idxs, dims_to_tensor_shape, dims_to_tensor_perm,
-    collapse_dims_super, collapse_dims_oper, einsum, Dimensions
+    collapse_dims_super, collapse_dims_oper, einsum, Dimensions,
+    to_tensor_rep, from_tensor_rep
 )
 
 
@@ -92,6 +93,10 @@ _Indices = collections.namedtuple('_Indices', ['base', 'permutation', 'shape'])
                           [2, 3, 0, 1, 6, 7, 4, 5],
                           (6, 8, 2, 4, 5, 7, 1, 3)),
                  id="super-oper"),
+    pytest.param(_Indices([[[1], [2, 3]], [[4, 5], [6]]],
+                          [2, 0, 1, 4, 5, 3],
+                          (2, 3, 1, 6, 4, 5)),
+                 id="super-oper-uneven"),
     pytest.param(_Indices([[[2, 4], [6, 8]], [1]],
                           [2, 3, 0, 1, 4], (6, 8, 2, 4, 1)),
                  id="operator-ket"),
@@ -114,6 +119,70 @@ class TestSuperOperatorDimsModification:
         assert dims_to_tensor_shape(indices.base) == indices.shape
 
 
+@pytest.mark.parametrize("dims", [
+    pytest.param([[2, 3], [1]], id="ket"),
+    pytest.param([[1, 1], [2, 3]], id="bra"),
+    pytest.param([[2, 3], [5, 7]], id="oper"),
+    pytest.param([[[2], [3]], [[5], [7]]], id="super-oper"),
+    pytest.param([[[2], [3]], [1]], id="oper-ket"),
+    pytest.param([[[1, 1], [1, 1]], [[2, 3], [5, 7]]], id="oper-bra"),
+    pytest.param([[[2], [3, 5]], [[7], [11]]], id="super-oper-uneven"),
+])
+class test_to_tensor_rep:
+    class mylist(list):
+        def at(self, i, j):
+            out = self.copy()
+            out[i] = j
+            return out
+
+    def _build_qobj(self, dims, idx=-1):
+        """
+        Create an Qobj with given dims which values only increment for the
+        `idx` dims.
+        """
+        if isinstance(dims[0], int):
+            ket = qutip.tensor([
+                (
+                    qutip.Qobj(np.arange(N))
+                    if i == idx else (qutip.Qobj(np.ones(N)))
+                ) for i, N in enumerate(dims)
+            ])
+            ket.dims = [ket.dims[0], [1]]
+            return ket
+
+        left = self._build_qobj(dims[0], idx)
+        right = self._build_qobj(dims[1], idx - len(flatten(dims[0])))
+        if left.isoper and left.type != 'scalar':
+            left = qutip.operator_to_vector(left)
+            left.dims = [left.dims[0], [1]]
+        if right.isoper and right.type != 'scalar':
+            right = qutip.operator_to_vector(right)
+            left.dims = [left.dims[0], [1]]
+        return left @ right.dag()
+
+    def test_tensor_shape(self, dims):
+        qobj = self._build_qobj(dims)
+        array = to_tensor_rep(qobj)
+        assert array.shape == tuple(flatten(dims))
+
+    def test_tensor_rep(self, dims):
+        N = len(flatten(dims))
+        for i in range(N):
+            qobj = self._build_qobj(dims, i)
+            array = to_tensor_rep(qobj)
+            slices = mylist([slice(j) for j in array.shape])
+            assert all(
+                np.all(array.__getitem__(*slices.at(i, j)) == j)
+                for j in range(array.shape[i])
+            )
+
+    def from_tensor_rep(self, dims):
+        qobj = self._build_qobj(dims)
+        array = to_tensor_rep(qobj)
+        back = from_tensor_rep(array, dims)
+        assert qobj == back
+
+
 class TestTypeFromDims:
     @pytest.mark.parametrize(["base", "expected"], [
         pytest.param([[2], [2]], 'oper'),
@@ -197,9 +266,9 @@ def test_dims_comparison():
     pytest.param("ij->ji", [qutip.sigmay()], qutip.sigmay().trans()),
     pytest.param("ij,ji", [qutip.sigmaz(), qutip.sigmaz()], 2),
     pytest.param("ijij", [qutip.tensor(qutip.thermal_dm(2,1), qutip.thermal_dm(2,1))], 1),
-    pytest.param("ikjl,jm->ikml", [qutip.tensor(qutip.sigmaz(), qutip.sigmaz()), 
+    pytest.param("ikjl,jm->ikml", [qutip.tensor(qutip.sigmaz(), qutip.sigmaz()),
                              qutip.sigmaz()], qutip.tensor(qutip.qeye(2), qutip.sigmaz())),
     pytest.param("ijkl->kjil", [qutip.tensor(qutip.sigmam(), qutip.sigmaz())], qutip.tensor(qutip.sigmap(), qutip.sigmaz()))
 ])
 def test_einsum(subscripts, operands, expected):
-    assert einsum(subscripts, *operands) == expected
\ No newline at end of file
+    assert einsum(subscripts, *operands) == expected