Also unexpectedly
np.exp(pm.Dirichlet.dist(a = np.array([2, 1])).logp(np.array([1, 0])).eval())
returns 0.0. But I assumed it to be
\frac{1}{B(a)} \prod_{i=1}^k x_i^{a_i-1} = \frac{1}{0.5} 1^1 \cdot 0^0 = 2.
The definition given here Multivariate — PyMC 5.10.0 documentation might has a typo since the last subscript k probably belongs to \theta like that:
f(x|\theta) = \frac{\Gamma( \sum_{i=1}^k \theta_i)}
{\prod_{i=1}^k \Gamma\left( \theta_i\right)}
\prod_{i=1}^{k-1} x_i^{\theta_i-1}
\left(1 - \sum_{i=1}^{k-1} x_i \right)^{\theta_k}
This would make sense because it states
Only the first k-1 elements of x are expected.
and since it should hold that \sum_{i=1}^k x_i = 1 the term in the brackets should equal x_k. But when tested
np.exp(pm.Dirichlet.dist(a = np.array([2, 1])).logp(np.array([0.5, 0.5])).eval())
returns 1.0 which is correct but
np.exp(pm.Dirichlet.dist(a = np.array([2, 1])).logp(np.array([.5, 0])).eval())
returns 0.
So I assume whenever there is a 0 in x it will return 0 even if the corespoding value in a is 1.