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$
2
  • $\begingroup$ "I was not successful" -- what went wrong? $\endgroup$
    – Chris K
    Mar 8, 2017 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, 2017 at 22:19

1 Answer 1

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$
1
  • $\begingroup$ Thank you very much for your suggestion @xzczd $\endgroup$
    – SAC
    Mar 10, 2017 at 16:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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