Needs["Graphics`PlotField`"];
l=100;
modell={u'[s]== u[s]-v[s]-u[s](u[s]^2+v[s]^2),
v'[s]== u[s]+v[s]-v[s](u[s]^2+v[s]^2),w'[s]== 0.03 w[s],u[0]==u0, v[0]==v0, w[0]==2};
p:=Module[{},loesung=NDSolve[modell,{u,v,w},{s,0,l}]; 
abb=ParametricPlot3D[Evaluate[{u[s],v[s],w[s]}/.loesung],{s,0,l},
DisplayFunction->Identity,
PlotPoints->500,PlotRange->{{-2,2},{-2,2},{0,1}}]
];
bild[1]=Table[{u0=z;v0=-2;p},{z,-2,2,2}];
bild[2]=Table[{u0=z;v0=2;p},{z,-2,2,2}];
bild[3]=Table[{u0=-2;v0=z;p},{z,-2,2,2}];
bild[4]=Table[{u0=2;v0=z;p},{z,-2,2,2}];
bild[5]=Table[{u0=z;v0=-0.1;p},{z,-0.5,0.5,0.5}];
bild[6]=Table[{u0=z;v0=0.1;p},{z,-0.5,0.5,0.5}];
bild[7]=Table[{u0=-0.1;v0=z;p},{z,-0.5,0.5,0.5}];
bild[8]=Table[{u0=0.1;v0=z;p},{z,-0.5,0.5,0.5}];
(*
feld=PlotVectorField[{z1-z2-z1(z1^2+z2^2), z1+z2-z2(z1^2+z2^2),0.01z3},{z1,-2,2,0.5},{z2,-2,2,0.5},
DisplayFunction->Identity];
*)
abb0=Show[Evaluate[Table[bild[i],{i,1,4}]],
DisplayFunction->$DisplayFunction,
PlotRange->{{-2,2},{-2,2},{2,5}},AxesLabel->{"x","y","z"},ImageSize->{400,300}];