Direction field for Newton’s second law:
({{{, ({{9.8-0.196v(t)=v'(t),  }, {v(0)=10,  }})}})

f[x_,y_]=9.8-0.196*y;
DSolve[{9.8-0.196*y[x]==y'[x],y[0]==10},y[x],x]
gr1=StreamPlot[{1,f[x,y]},{x,0,150},{y,0,150},Frame->False,Axes->True,AspectRatio->1/GoldenRatio,AxesLabel->{t,v[t]}];
gr2=Plot[50.` (E^(-0.196` x)) (-0.8`+1.` (E^(0.196` x))),{x,-25,150},PlotStyle->{Orange,Dashed,Thick}];
Show[gr1,gr2]

Direction field for equation
x'(t)=2x(t)

f[x_,y_]=2*y;
gr1=StreamPlot[{1,f[x,y]},{x,-2,2},{y,-2,2},Frame->False,Axes->True,AspectRatio->1/GoldenRatio,AxesLabel->{t,x[t]}];
gr2=Plot[-3*(E^(2 x)),{x,-2,2},PlotStyle->{Orange,Dashed,Thick}];
Show[gr1,gr2]

Nonautonomous problem
x’(t)=2x(t)t

f[x_,y_]=x*y;
DSolve[{y'[x]==y[x]*x,y[0]==1},y[x],x]
gr1=StreamPlot[{1,f[x,y]},{x,-2,2},{y,-2,2},Frame->False,Axes->True,AspectRatio->1/GoldenRatio];
gr2=Plot[(E^((x^2)/2)),{x,-2,2},PlotStyle->{Orange,Dashed,Thick}];
Show[gr1,gr2]

An equation with a nonunique solution

f[x_,y_]=Sqrt[y];
DSolve[{y'[x]==Sqrt[y[x]],y[0]==0},y[x],x]
gr1=StreamPlot[{1,f[x,y]},{x,-2,2},{y,0,2},Frame->False,Axes->True,AspectRatio->1/GoldenRatio];
gr2=Plot[(x^2)/4,{x,-2,2},PlotStyle->{Orange,Dashed,Thick}];
gr3=Plot[0,{x,-5,3},PlotStyle->{Green,Dashed,Thick}];
Show[gr1,gr2,gr3]

An equation with a steady state

Clear[x,y,f];
f[x_,y_]=2-y;
y0val=2;
FullSimplify[DSolve[{y'[x]==2-y[x],y[0]==y0},y[x],x]]
gr1=StreamPlot[{1,f[x,y]},{x,-2,2},{y,-2,3},Frame->False,Axes->True,AspectRatio->1/GoldenRatio];
gr2=Plot[(E^(-x)) (-2+2 (E^x)+y0val),{x,-2,2},PlotStyle->{Orange,Dashed,Thick}];
Show[gr1,gr2]

Another equation with a steady state
x’(t)=6-2*x(t)
x(0)=x0

Clear[x,y,f];
f[x_,y_]=6-2*y;
y0val=10;
DSolve[{y'[x]==6-2*y[x],y[0]==y0},y[x],x]
gr1=StreamPlot[{1,f[x,y]},{x,-2,2},{y,-2,10},Frame->False,Axes->True,AspectRatio->1/GoldenRatio];
gr2=Plot[(E^(-x)) (-2+2 (E^x)+y0val),{x,-2,2},PlotStyle->{Orange,Dashed,Thick}];
Show[gr1,gr2]

A phase diagram for a complicated system of differential equations

splot=StreamPlot[{y,(1-x^2) y-x},{x,-4,4},{y,-3,3},StreamColorFunction->"Rainbow"];
Manipulate[Show[splot,ParametricPlot[Evaluate[First[{x[t],y[t]}/.NDSolve[{x'[t]==y[t],y'[t]==y[t] (1-x[t]^2)-x[t]+0.5 Cos[1.1 t],Thread[{x[0],y[0]}==point]},{x,y},{t,0,T}]]],{t,0,T},PlotStyle->Red]],{{T,20},1,100},{{point,{3,0}},Locator},SaveDefinitions->True]

Phase diagrams for different versions of Romeo and Julliette problems.

(*m={{1,0},{0,1}}*)
m={{2,0},{0,-1}};
m//MatrixForm

StreamPlot[Evaluate[m.{x,y}],{x,-1,1},{y,-1,1}, StreamScale->Large, StreamColorFunction->"Rainbow",FrameLabel->{R,J},ImageSize->Small]

m.{x,y}

m={{0,1},{-1,0}};
m//MatrixForm
splot=StreamPlot[m.{x,y},{x,-4,4},{y,-3,3},FrameLabel->{R,J},StreamColorFunction->"Rainbow"];
Manipulate[Show[splot,ParametricPlot[Evaluate[First[{x[t],y[t]}/.NDSolve[{m.{x[t],y[t]}=={x'[t],y'[t]},Thread[{x[0],y[0]}==point]},{x,y},{t,0,T}]]],{t,0,T},AxesLabel->{R,J},PlotStyle->Red]],{{T,20},1,100},{{point,{3,0}},Locator},SaveDefinitions->True]

m.{x,y}

m={{-1,2},{1,-1}};
m//MatrixForm
splot=StreamPlot[m.{x,y},{x,-4,4},{y,-3,3},StreamColorFunction->"Rainbow"];
eigen1=ParametricPlot[{-(Sqrt[2])*u,1*u},{u,-3,3}];
eigen2=ParametricPlot[{(Sqrt[2])*u,1*u},{u,-3,3}];
Manipulate[Show[eigen1,eigen2,splot,ParametricPlot[Evaluate[First[{x[t],y[t]}/.NDSolve[{m.{x[t],y[t]}=={x'[t],y'[t]},Thread[{x[0],y[0]}==point]},{x,y},{t,0,T}]]],{t,0,T},PlotStyle->Red]],{{T,20},1,100},{{point,{3,0}},Locator},SaveDefinitions->True]
Eigenvalues[m]
Eigenvectors[m]

m={{-2,1},{1,-2}};
m//MatrixForm
splot=StreamPlot[m.{x,y},{x,-4,4},{y,-3,3},FrameLabel->{R,J},StreamColorFunction->"Rainbow"];
Manipulate[Show[splot,ParametricPlot[Evaluate[First[{x[t],y[t]}/.NDSolve[{m.{x[t],y[t]}=={x'[t],y'[t]},Thread[{x[0],y[0]}==point]},{x,y},{t,0,T}]]],{t,0,T},AxesLabel->{R,J},PlotStyle->Red]],{{T,20},1,100},{{point,{3,0}},Locator},SaveDefinitions->True]