The common way is to use a link function that constrains your expectation to be positive, such as exp or softplus.
mu = pm.math.exp(fn(a, b, c, x1, x2))
Alternatively, you can add an arbitrary penalty term with Potential, but NUTS will usually struggle with hard boundaries like that (so not recommended if first approach is fine)
mu = fn(a, b, c, x1, x2)
pm.Potential("negative_penalty", pm.math.switch(mu < 0, -np.inf, 0))