code: def density_matrix(self):
          
          size = 2 ** len(self._qubit_map)
          return np.reshape(self._density_matrix, (size, size))
docs: Returns the density matrix at this step in the simulation.
  
          The density matrix that is stored in this result is returned in the
          computational basis with these basis states defined by the qubit_map.
          In particular the value in the qubit_map is the index of the qubit,
          and these are translated into binary vectors where the last qubit is
          the 1s bit of the index, the second-to-last is the 2s bit of the index,
          and so forth (i.e. big endian ordering). The density matrix is a
          `2 ** num_qubits` square matrix, with rows and columns ordered by
          the computational basis as just described.
  
          Example:
               qubit_map: {QubitA: 0, QubitB: 1, QubitC: 2}
               Then the returned density matrix will have (row and column) indices
               mapped to qubit basis states like the following table
  
                      | QubitA | QubitB | QubitC
                  :-: | :----: | :----: | :----:
                   0  |   0    |   0    |   0
                   1  |   0    |   0    |   1
                   2  |   0    |   1    |   0
                   3  |   0    |   1    |   1
                   4  |   1    |   0    |   0
                   5  |   1    |   0    |   1
                   6  |   1    |   1    |   0
                   7  |   1    |   1    |   1