Submanifolds of constant entanglement entropy In a Hilbert space $A \otimes B$, a density matrix $\rho: A \rightarrow A$ has associated with it an entanglement entropy $S(\rho)$.
Question: What is the description of the collection of all the states $\rho'$ such that $S(\rho) = S(\rho')$?

In case the answer isn't obvious, here are some sub-questions.
Notation


*

*Let $H(A)$ denote the space of density matrices $A \to A$. This is a topological space via embedding into the space of linear automorphisms.

*For $\rho \in H(A)$, denote by $[\rho] \subset H(A)$ equivalence classes of the relation given by equivalence of entanglement entropy.
$$\rho \sim \rho' \iff S(\rho) = S(\rho')$$
Sub-questions


*

*Is there a transitive group action on $[\rho]$? (I think I read somewhere unitary maps preserve entanglement entropy. Are they transitive in these equivalence classes?)

*Is the quotient map $H(A) \xrightarrow{\pi} (H(A)/\sim)$ a fibration?

*Topologically, are $[\rho]$s compact? Do they have boundary?

*Is $[\rho]$ connected? What does a connected component look like?

*Are there infinitesmall deformations of a state $\rho$ that preserve entropy? Equivalently, the classes $[\rho]$ probably inherit a smooth structure, what are the tangent spaces? In the same vein, is there an obstruction theory for lifting these 1st order deformations to higher order deformations?

*The following is non-rigorous. One often computes entanglement entropy in a continuum lattice scenario, e.g. in CFTs. Here one considers, say, a compact orientable surface, choose a closed curve $\gamma$ in the trivial homology class, and consider entanglement entropy of a state in "the space of states bounded by the curve". There are standard CFT techniques (due to Cardy probably?) for computing these, e.g. the replica trick. In these scenarios, in addition to the questions above, one can ask the following. Fixing a state, are there deformations of $\gamma$ that preserve entanglement entropy?

*More generally, in other scenarios one can compute some notion of entropy (e.g., Shannon/Von Neuman entropy in classical probability), are there answers to the above questions?

*Speculative: if answers to all the above are positive, one could further ask for a symplectic structure on $H(A)$ under which the quotient is a Lagrangian fibration with Lagrangian fibers. Is anything like this known?

Exercise/Intuition
Bloch Sphere (edit: this was done below by @Noiralef):
It should be straightforward to compute these level surfaces of the entanglement entropy function on the Bloch sphere, which could provide intuition for the more general situation. Though it was mentioned in the comments that in this special case one does expect the unitaries to be transitive, which won't be the case in general.
Shannon Entropy
There's another nice computation one can do. Consider a classical probabilistic model on a finite number of variables, the state space in $n$-variables is an $n-1$ simplex in Euclidean space. Observe the Shannon entropy has a $S_{n-1}$ symmetry by swapping any pair of variables. Hence when performing a quotient of the type described above one is bound to end up with an orbifold (or a DM-stack).
Actually looks like these in this case are Morse functions too, so as usual the fibers of the map degenerate at critical points. Together with the observation above looks like we should be doing equivariant Morse theory.
 A: This will only be a partial answer, because I don't even know all of those words. But it seems to me that you are "thinking too complicated" about the topic, so I hope this will be helpful.
First of all, you are talking about "entanglement entropy". As was already mentioned in the comments, that is an overly fancy word here: Your space $B$ plays no role whatsoever, we are simply talking about the Hilbert space $A$ and the (von Neumann) entropy $S: H(A) \to \mathbb R$.

Let us first talk about a qubit, $A = \mathbb C^2$.
$H(A)$ can be parametrized with the Bloch vector, i.e., we write a density matrix $\rho \in H(A)$ as
$$ \rho = \frac 1 2 \begin{pmatrix} 1 + z & x - \mathrm i y \\ x + \mathrm i y & 1 - z \end{pmatrix} . $$
The Bloch vector $\vec r = (x, y, z)^t$ lies in the Bloch sphere $B^2 = \{ \vec r \in \mathbb R^3: | \vec r | \leq 1 \}$.
Some simple algebra tells us that
$$ S(\rho) = -\operatorname{tr}(\rho \log \rho) = -\frac 1 2 \log \frac{1 - |\vec r|^2}{4} - \frac{|\vec r|}{2} \log \frac{1 + |\vec r|}{1 - |\vec r|} . \tag{1} $$
$S(\rho)$ only depends on the length of the Bloch vector $\vec r$, and you can easily check that (1) is a strictly decreasing function of $|\vec r|$ (which is $\log 2$ at the maximally mixed state $\vec r = 0$ and $0$ at pure states with $|\vec r| = 1$).
Therefore an equivalence class $[\rho]$ is exactly a shell $S^2(r) = \{ \vec r \in \mathbb R^3: |\vec r| = r \}$ (with $0 \leq r \leq 1$).
They are compact, without boundary, connected, and $\operatorname{SO}(3)$ acts transitively on them.

We continue with the general case $A = \mathbb C^N$.
$H(A)$ can be very complicated geometrically:
If we parametrize the density matrix with a generalized Bloch vector, the space of valid parameters will become much more complicated.
Topologically, there is not much difference though.
$H(A)$ is still a convex compact space, with the boundary consisting of the pure states.
Let
$$ H_S(A) = \{ \rho \in H(A): S(\rho) \geq S \} , $$
then $H(A) = H_0(A) \supset \cdots \supset H_{\log N}(A)$ is a nested family of convex compact spaces (because $S$ is a concave function), their boundaries are your equivalence classes.
All shells / equivalence classes $\partial H_S(A)$ -- except the innermost $H_{\log N}(A) = \{ \frac 1 N \mathbb 1 \}$, which is only the maximally mixed state -- are homeomorphic to the outermost shell $\partial H(A)$ of pure states. To construct the homeomorphism note that the ray from $\frac 1 N \mathbb 1$ through $\rho \in \partial H_S(A)$ hits $\partial H(A)$ exactly once.
So, all equivalence classes are the same as the space of pure states $\partial H(A)$.
A pure state is given by a vector $0 \neq \psi \in A$ up to equivalence $\psi \sim \lambda \psi\,$ ($0 \neq \lambda \in \mathbb C$), therefore $\partial H(A)$ is just projective $\mathbb C P^{N-1}$.
