I have the following system of equations
$$ \rho V(u, \epsilon^i) = F(u, \epsilon^i) + V_u(u, \epsilon^i)g(u, V(u, \epsilon^i) + \lambda^i \left(V(u, \epsilon^{-i}) - V(u, \epsilon^i)\right)$$
with $i \in \{0,1\}$, so that $\epsilon$ can take two values and $g(u, V)$ known and initial $u(0), \epsilon$ given. $g(u,V)$ is the law of motion for $u$, $\dot u(u, \epsilon^i)$ if you will - but it depends on the value of $V$. This looks like a HJB, but crucially, there is no optimization involved.
This is a macro-labor model, where the amount of job openings depends on the value of a filled position ($V$), which is why the law of motion for $u$ (unemployment) depends on the equilibrium value of $V$.
How could I compute (numerically?) the value of $V$ that solves this ODE?
