I would like to build a system of differential equations in $\mathbb{R}^2$ that has a set of pre-determined invariant curves. For instance, given the three following functions \begin{equation} g_1(x.y)=x\qquad g_2(x,y)=y\qquad g_3(x,y)=x-1 \end{equation} I'd like the system to have $g_1=0,g_2=0$ and $g_3=0$ as invariant curves (lines, in this case).
If the system is given by \begin{equation} (\dot{x},\dot{y})=F(x,y) \end{equation} I thought that the condition to set had to be \begin{equation} \nabla g_i(x,y)\cdot F(x,y)=0\qquad i=1,2,3 \end{equation} but I'm not sure as these conditions yield unreasonable constraints...any help?
If you have vector function $\boldsymbol f$ and scalar function $g$, then $$ \boldsymbol F = (\boldsymbol \nabla g)^2\boldsymbol f-(\boldsymbol f\boldsymbol \nabla g)\boldsymbol\nabla g $$ will be always collinear with vector field $\boldsymbol\nabla g$.
This gives us a hint how to build what you want. For some arbitrary functions $\boldsymbol f_i$: $$ F = \left(\boldsymbol f_0+\sum_{i=1}^{3}\frac{ (\boldsymbol\nabla g_i)^2\boldsymbol f_i-(\boldsymbol f_i\boldsymbol\nabla g_i)\boldsymbol\nabla g_i }{g_i}\right)g_1g_2g_3. $$
For example, for $f_0=f_1=f_2=f_3=(x^2-y^2,x+y)$ and your $g_i$: $$ F=\left(x(x-1) (y+1) \left(x^2-y^2\right),\left(x^2+x-1\right) y (x+y)\right) $$