Consider the following system defined on the open interval (-1, 1):
$y_1' = c y_3 \\ y_2' = c y_4 \\ y_3' = -f(y_1, y_2)y_2 \\ y_4' = f(y_1, y_2)y_1 $
given $ y_3(-1) = 0 = y_3(1) \\ y_4(-1) = 1, y_4(1) = -1 $
Note that the system cannot be reduced to the standard form of two point BVP in $y_3, y_4$ by eliminating $y_1, y_2$.
What methods/packages can I use to solve this system ?