Basically, I have a set of differential equations that I need to solve for exactly 100 different initial conditions (given as lists for each initial condition), and then plot each solution.

Here is some sample code where I have set vrad, vtan, and deltaR (arrays of initial conditions) to an array of length two. So, given the arrays vrad, vtan, deltaR (our initial conditions) I want to be able to essentially do what this code does but for the array of solutions. Cheers!

Edit: I think I've nearly done it, I just need Table to not iterate through every tuple, but instead by index, anyone know how to do this?

(* Scaling Quantities *)
V = 200;
R = 10^4;
(* Random Quantities *)
vrad = {0, 5};
vtan = {0, 5};
deltaR = {0, 5};
(* Converting to dimensionless quantities *)
vRadial = (V + vrad)/V;
vTangential = (V + vtan)/V;
r0 = (10^4 + deltaR)/R;
L = r0*vTangential;
(* numerical solution *)
s = Partition[
  Flatten@Table[
    NDSolve[{r''[t] == r[t]*ϕ'[t]^2 - 1/r[t], ϕ'[t] == d/
       r[t]^2, ϕ[0] == a, r[0] == b, 
      r'[0] == c}, {r, ϕ}, {t, 0, 200}], {a, vTangential/r0}, {b,
      r0}, {c, vRadial}, {d, L}], 2]
(* Plotting the solution *)
ParametricPlot[
 Evaluate[{r[t]*Cos[ϕ[t]], r[t]*Sin[ϕ[t]]} /. s], {t, 0, 
  2*Pi}, GridLines -> Automatic, Frame -> True]