System von nicht linearen DGL nummerisch berechnen

Hallo,

ich habe kürzlich in c++ lotka volterra Beute/Jäger
System geschrieben. Dabei habe ich RungeKutta 4 Ordnung
zum Berechnen des gekoppelten Systems benutzt.
Nun würde ich gerne in Matlab oder Octave dasselbe
modellieren. Später vielleicht auch konkurriernde
Populationen, was im Prinzip auch auf gekoppeltes
System von DGLs hinausläuft.
Zustandsraummodell mit A,B,C,D Matrizen scheiden leider aus.
x1 und x2 sind nicht linear, sondern durch Multiplikation verbunden,
damit kann A nicht formuliert werden.
Wie kann ich es also lösen?

Grüsse, Daniel