I am trying to set up a recursive function that generates n number of differential equations for Subscript[y, n][t]
This function almost works.
Table[{Subscript[y0, j] = 1}, {j, 50}];(*initial conditions for Subscript[y, n] assuming n<=50*) vars[n_] := {x, Table[Subscript[y, j], {j, n}]}; sol[0][T_, b_, d_, r_, n_] := sol[0][T, b, d, r, n] = Flatten[vars[n]] /. NDSolve[Flatten@{Join[{x'[t] == 10^7 r - d x[t] - b*x[t]*(Sum[Subscript[y, k][t], {k, n}])}, Table[Subscript[y, j]'[t] == -d Subscript[y, j][t] + b x[t]*Subscript[y, j][t], {j, n}], Flatten[Join[{x[0] == 0}, Flatten[Table[{Subscript[y, j][0] == Subscript[y0, j]}, {j, n}]]]]]}, Flatten[vars[n]], {t, 0, T}][[1]] sol[i][T_, b_, d_, r_, n_] := sol[i][T, b, d, r, n] = Flatten[vars[n]] /. NDSolve[Flatten@{Join[{x'[t] == 10^7 r - d x[t] - b*x[t]*(Sum[Subscript[y, k][t], {k, n}])}, Table[Subscript[y, j]'[t] == -d Subscript[y, j][t] + b x[t]*Subscript[y, j][t], {j, n}], Flatten[Join[{x[0] == 0},(*next bit seems to be the problem*) Flatten[Table[{Subscript[y, j][0] == sol[i - 1][T, b, d, r, n][[j + 1]][T]}, {j, n}]]]]]}, Flatten[vars[n]], {t, 0, T}][[1]] The initial condition sol[0][T, b, d, r, n] works as expected and returns the interpolating functions:
T = 4; b = 10^-7; d = 0.25; r = 0.2; n = 4; sol[0][T, b, d, r, n] And also returns the value at t = T e.g. sol[0][T, b, d, r, n][[2]][4]
sol[i][T, b, d, r, n] does not work for i>0.
It seems that the problem is where the solution from the previous iteration i-1 is used to set the initial conditions for current iterate i, as marked in the code.
I imagine this will be trivial to troubleshoot for someone on here. Any advice is much appreciated.
