Issue with ParametricNDSolve and Table

Question

I have a differential equation, I plot it for different initial values by using ParametricNDSOlve.

Code - answer = ParametricNDSolve[{x'[t]==-3x[t] - y[t],y'[t] == x[t],x[0]==a,y[0]==a},{x,y},{t,0,100},{a}].

Plotting - pp = ParametricPlot[Evaluate@Table[{x[a][t],y[a][t]}/.answer,{a,-2,2,0.4}],{t,0,40},PlotRange->All,PlotLegends->Range[-2,2,0.4]]; pp/.Line[x_]:>{Arrowheads[Table[.04, {4}]], Arrow[x]}

Issue with above code that all initial conditions are of the form (a,a), I need to solve for (a,b). So I use two parameters for in the ParametricNDSolve - answer = ParametricNDSolve[{x'[t]==-3x[t] - y[t],y'[t] == x[t],x[0]==a,y[0]==b},{x,y},{t,0,100},{a,b}]. Till solving part is fine. Issue comes while plotting as I use the command Table while generating multiple values of initial position at time=0. If it was a single variable Table was generating a 1-D array for the ParametricPlot to use. Now When I use something like ParametricPlot[Evaluate@Table[{x[a,b][t],y[a,b][t]}/.answer, {a,-2,2,0.5},{b-2,2,0.5}],{t,0,40}], Table gives a matrix for parameters(what like seems to me) but what I need is something like a 2-D matrix like one value of a and one value of b. Something like (-2,-2),(-2,-1.5),(-2,1)...(-2,2),(-1.5,-2),(-1.5,-1.5) ....

Also like If someone could explain what is happening pp/.Line[x_]:>{Arrowheads[Table[.04, {4}]], Arrow[x]} would also be very nice.

Answer 1

Clear["Global`*"]

answer1 = ParametricNDSolve[
   {x'[t] == -3 x[t] - y[t], y'[t] == x[t], x[0] == a, y[0] == a},
   {x, y}, {t, 0, 100}, {a}];

pp1 = ParametricPlot[
   Evaluate[{x[a][t], y[a][t]} /. answer1],
   {a, -2, 2}, {t, 0, 40},
   PlotRange -> All,
   PlotPoints -> 25,
   MaxRecursion -> 3,
   BoundaryStyle -> None];

EDIT: Corrected figure

pp2 = ParametricPlot[Evaluate@
    Table[Tooltip[{x[a][t], y[a][t]}, a] /. answer1,
     {a, 2, -2, -0.4}], {t, 0, 40},
   PlotRange -> All,
   PlotLegends -> LineLegend[Range[2, -2, -0.4],
     LegendLabel -> Style["a", 14, Bold],
     LegendMarkerSize -> 20,
     LegendLayout -> {"Column", 1}]];

Show[pp1, pp2 /.
  Line[x_] :> {Arrowheads[Table[.04, {4}]], Arrow[x]}]

answer2 = ParametricNDSolve[
   {x'[t] == -3 x[t] - y[t], y'[t] == x[t], x[0] == a, y[0] == b},
   {x, y}, {t, 0, 100}, {a, b}];

pp3 = ParametricPlot[Evaluate[
    Table[{x[a, b][t], y[a, b][t]}, {b, -2, 2, 0.5}] /. answer2],
   {a, -2, 2}, {t, 0, 40},
   PlotRange -> All,
   PlotPoints -> 25,
   MaxRecursion -> 3,
   BoundaryStyle -> None];

pp4 = ParametricPlot[Evaluate@Flatten[
     Table[
      Tooltip[{x[a, b][t], y[a, b][t]}, {a, b}] /. answer2,
      {b, 2, -2, -0.5}, {a, 2, -2, -0.5}], 1], {t, 0, 40},
   PlotStyle -> (ColorData[97] /@ Range[9]),
   PlotRange -> All,
   PlotLegends -> {LineLegend[Range[2, -2, -0.5],
      LegendMarkerSize -> 20,
      LegendLayout -> {"Column", 1},
      LegendLabel -> Style["a", 14, Bold]],
     SwatchLegend[Range[-2, 2, 0.5],
      LegendMarkerSize -> 10,
      LegendLayout -> {"Row", 1},
      LegendLabel -> Style["b", 14, Bold]]}];

Show[pp3, 
 pp4 /.
  Line[x_] :> {Arrowheads[Table[.04, {4}]], Arrow[x]}]