known change of variables leads to the solution. Indeed, and divide by $f_i^{1-\gamma_i}$ (after noting that $\gamma_i > 0$) and use the substitution $u_i(t) = f_i^{\gamma_i}(t)$, so that (\[eqf\]) becomes the linear differential equation $$\frac{1}{\gamma_i} u_i' + a_i u_i + b_i = 0,$$ with terminal condition $u_i(T) = 1$. This linear equation admits the unique solution $$u_i(t) = e^{\gamma_i \int_t^T a_i(s)ds} + \int_t^T \gamma_i b_i(s) e^{- \gamma_i \int_t^s a_i(r)dr}ds.$$ Note that $b_i$ is positive everywhere, and thus so is $u_i$. Hence, $u_i^{1/\gamma_i}$ is well defined, and the unique solution to (\[eqf\]) is given by $$\label{solBern}
f_i(t) = \left(e^{\gamma_i \int_t^T a_i(s)ds} + \int_t^T \gamma_i b_i(s) e^{- \gamma_i \int_t^s a_i(r)dr}ds \right)^{1/\gamma_i}.$$ Substituting this solution into the ansatz yields the solution $v(x,y,t)$ of the HJB equation, as long as we check that $v_{xx} < 0$ and $v_x > 0$. But this is straightforward: $$\begin{aligned}
v_x(x,y,t) &= \epsilon_i (1 - \theta_i/n) x^{- 1 / \gamma_i} y^{- \theta_i( 1 - 1/\delta_i)} f_i(t) > 0, \\ 
v_{xx}(x,y,t) &= - \frac{\epsilon_i}{\gamma_i} (1 - \theta_i/n) x^{- 1 / \gamma_i - 1} y^{- \theta_i( 1 - 1/\delta_i)} f_i(t) < 0,\end{aligned}$$ where we again use $\gamma_i > 0$. Therefore, in terms of $f_i$, we may express the optimal controls from 