0
$\begingroup$

I have an ODE system which solves of n variables, with initial conditions defined using the previous differential equation solution of n-1 variables and with an initial condition for the last variable (which randomly depends on one of the previous variables r[n_]). In addition, the system has a other random variable (\[Xi][j_, k_] different of r[n_]). The code of the system is shown below:

nmax = 5;(*max number of equations*)
T = 70;(*Max Time*)
\[Xi][j_, k_] := RandomReal[];(*random variable*)
r[n_] := RandomInteger[{1, 
   n}];(*set up unknown vars and differential equations,for n \
species*)
vars := Table[Subscript[x, j], {j, n}];
eqns := Table[
   Subscript[x, j]'[t] == 
    Subscript[x, j][
      t] (1 - Subscript[x, j][
        t] - (Sum[\[Xi][j, k] Subscript[x, k][t] Boole[k != j], {k, 
          n}])), {j, n}];
ics = {Subscript[x, 1][0] == 0.7};(*initial ICs*)
(*main loop*)
Do[(*solve for n species*)
 sol[n] = NDSolve[{eqns, ics}, vars, {t, 0, T}][[1]];
 (*set up ICs for n=n+1 species*)
 ics = Join[
   Table[Subscript[x, j][0] == 
     Evaluate[Subscript[x, j][T] /. sol[n]], {j, 
     n}], {Subscript[x, n + 1][0] == 
     Evaluate[0.01 Subscript[x, r[n]][T] /. sol[n]]}];, {n, nmax}]

I try to do two things:

i) Do 1000 different process as described above;

I thought about adding a Do as:

Do[Do[(*solve for n species*)
 sol[n] = NDSolve[{eqns, ics}, vars, {t, 0, T}][[1]];
 (*set up ICs for n=n+1 species*)
 ics = Join[
   Table[Subscript[x, j][0] == 
     Evaluate[Subscript[x, j][T] /. sol[n]], {j, 
     n}], {Subscript[x, n + 1][0] == 
     Evaluate[0.01 Subscript[x, r[n]][T] /. sol[n]]}];, {n, nmax}],{pp,1000}]

But I was not successful.

ii) Obtain a list with all the values of the variable \[Xi][j, k] drawn in the 1000 operations.

Can someone help me?

$\endgroup$
  • $\begingroup$ "I was not successful" -- what went wrong? $\endgroup$ – Chris K Mar 8 '17 at 20:37
  • $\begingroup$ @ChrisK The errors: NDSolve::ndode: Input is not an ordinary differential equation. ReplaceAll::reps: {(Subscript[x, 1]^[Prime])[t]==(1. -Subscript[x, 1][t]) Subscript[x, 1][t]} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing. ReplaceAll::reps: {Subscript[x, 1][0]==0.223725,Subscript[x, 2][0]==1.39449*10^-8,Subscript[x, 3][0]==0.000173822,Subscript[x, 4][0]==0.715508,Subscript[x, 5][0]==0.810497,Subscript[x, 6][0]==0.00223725} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing. $\endgroup$ – SAC Mar 8 '17 at 22:19
3
$\begingroup$
  1. Your attempt in (i) fails because you forgot to initialize ics for every round of repeation.
  2. Generating the list of ξ[j, k] before the summation will make the collection of it easier.
  3. It's better to define vars and eqns as functions of n.

The following is my solution:

nmax = 5; T = 70; repeat = 1000;
r[n_] := RandomInteger[{1, n}];
vars[n_] := Table[x[j], {j, n}]; 
eqns[n_, ξ_] := 
 Table[x[j]'[t] == 
   x[j][t] (1 - x[j][t] - Sum[ξ[[j, k]] x[k][t] Boole[k != j], {k, n}]), {j, n}];

ans = Table[Module[{ics = {x[1][0] == 0.7}},
    Table[
     Module[{xi = RandomReal[1, {n, n}], sol, icvalue},
      sol = NDSolveValue[{eqns[n, xi], ics}, vars[n], {t, 0, T}];
      icvalue = Through@sol[T];
      ics = {Thread[Table[x[j][0], {j, n}] == icvalue], 
        x[n + 1][0] == 0.01 icvalue[[r@n]]};
      {sol, xi}], {n, nmax}]],
   {repeat}];
$\endgroup$
  • $\begingroup$ Thank you very much for your suggestion @xzczd $\endgroup$ – SAC Mar 10 '17 at 16:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.