How can one integrate over the unit cube, subject to certain (quantum-information-theoretic) constraints? To begin, we have two constraints
\begin{equation}
C1=x>0\land z>0\land y>0\land x+2 y+3 z<1
\end{equation}
and
\begin{equation}
C2=x>0\land y>0\land x+2 y+3 z<1\land x^2+x (3 z-2 y)+(y+3 z)^2<3 z.
\end{equation}
$C1$ ensures the nonnegative-definiteness of a class of $9 \times 9$ ("two-qutrit") density matrices ($\rho$).
$C2$ also ensures this, as well as the nonnegative-definiteness of the "partial transpose" ($\rho^{PT}$) of $\rho$.
The integration--subject to $C1$--of the value 36 over the unit cube $\{x,y,z\} \in [0,1]^3$ yields 1.
The integration--subject to $C2$--of the value 36 over the unit cube yields (the Hilbert-Schmidt positive partial transpose probability) 
$\frac{8 \pi}{27 \sqrt{3}} \approx 0.537422$.
Now, we are interested in similarly enforcing both $C1 \land C3$ (yielding an "entanglement probability") and $C2 \land C3$ (yielding a "bound-entanglement probability"),
where, the entanglement constraint $C3$ is 
\begin{equation}
b \left(-x \left(a^2-a (b+2)+2 b+1\right)+2 y (a (-a+b+2)+b-1)-3 z (a-b-1)
   (a+b-1)+(a-1)^2\right)<0,
\end{equation}
with its three parameters $a,b,c$ subject to
\begin{equation}
C4= b>0\land c>0\land 0<a<1\land a+b+c=2\land (a-1)^2=b c.
\end{equation}
Now, I suspect these last two problems are too difficult to resolve in their full generality (leaving $a,b,c$ unspecified).
But, to begin, if we take
\begin{equation}
\{a,b,c\}=\left\{\frac{1}{4} \left(3-\sqrt{5}\right),\frac{1}{2},\frac{1}{4}
   \left(3+\sqrt{5}\right)\right\},
\end{equation}
the integration of 36 over the unit cube subject to $C1 \land C3$ yields 
$\frac{5}{132} \left(5+\sqrt{5}\right) \approx 0.274093$.
Alternatively, for 
\begin{equation}
\{a,b,c\} = \left\{\frac{1}{3},\frac{1}{3},\frac{4}{3}\right\},
\end{equation}
the integration of 36 over the unit cube subject to $C1 \land C3$ yields 
$\frac{125}{486}  \approx 0.257202$.
However, I have not been so far able to obtain the counterparts for these last two results for $C2 \land C3$. (Using numerical integration, we get the much lower values of 0.001497721920258410 and 0.003256122941383665, respectively.)
A set of values of $a,b,c$ which satisfy $C4$ are 
\begin{equation}
\left\{\frac{2}{3} (\cos (\alpha )+1),\frac{2}{3} \left(-\frac{1}{2} \sqrt{3} \sin
   (\alpha )-\frac{\cos (\alpha )}{2}+1\right),\frac{2}{3} \left(\frac{1}{2} \sqrt{3}
   \sin (\alpha )-\frac{\cos (\alpha )}{2}+1\right)\right\},
\end{equation}
for $\frac{\pi}{3} \leq \alpha \leq \frac{5 \pi}{3}$.
So, I would like to obtain results of integration of the value 36 over $[0,1]^3$ of $C1 \land C3$ and $C2 \land C3$ for either specific values of 
$a,b,c$, satisfying $C4$, or even without particular values being specified.
 A: Well, maybe the second comment of LSpice made me fully realize that the constraint $C4$ means that there is only one degree of freedom between $a,b,c$. So, as a start I took $a=\frac{1}{3}$. Then, using the Mathematica GenericCylindricalDecomposition command and choosing the order of integration of $x,y,z$ to select the simplest (as measured by LeafCount) output, I was able to obtain an entanglement probability of 
$\frac{125}{486}=\frac{5^3}{2 \cdot 3^5} \approx 0.257202$.
Additionally, using the same approach, the bound-entanglement probability proved to be 
\begin{equation}
\frac{-204+56 \sqrt{3} \pi +7 \log (7)-336 \sqrt{3} \csc ^{-1}\left(2
   \sqrt{7}\right)}{1134} \approx 0.00325612.
\end{equation}
So, I'll now try to use other values of $a$--rather than $\frac{1}{3}$ (and extend this answer, if I so succeed, and the results are of interest).
I thought that the problem of getting results for general $0 \leq a \leq 1$ would be much more formidable, as the resulting forms that $C1 \land C3$ and $C2 \land C3$ take seemed quite involved.
However, much to my surprise, the calculation for $C1 \land C3$ for the entanglement probability as a function of $a$ greatly simplified, yielding
\begin{equation}
-\frac{(a-2)^3}{9 a^2-30 a+27}.
\end{equation}
For $a=\frac{1}{3}$, the function gives the above-reported $\frac{125}{486}$.
Further, Nicholas Tessore was able to find the formula for $C2 \land C3$, giving the bound-entanglement probability.
It took the form
\begin{equation} \label{Tessore}
-\frac{A+B}{54 (4-3 a)^{3/2} (2 a-3)} ,   
\end{equation} 
where
\begin{equation}
A=8 \sqrt{12-9 a} \left(6 a^2-17 a+12\right) \cos ^{-1}\left(\frac{a (3 a-8)+6}{6-4
   a}\right)
\end{equation}
and
\begin{equation}
B=3 \sqrt{a} \left(2 \left(9 a^3-57 a^2+108 a-64\right)+3 (3-2 a) a \log (9-6 a)\right).
\end{equation}
So, a complete answer has now been successfully given to the question put.
