0
$\begingroup$

I have an ODE system that increases in size according to the rules

n = 5;
T = 50;
nu = 0.05;
vars = Table[Subscript[x, j][t], {i, n}, {j, i}];
eqns = Table[{Subscript[x, j]'[t] == 
    Subscript[x, j][
      t] (1 - Subscript[x, j][t] - 
       nu (Sum[ Subscript[x, k][t] Boole[k != j], {k, i}]) ), 
   Subscript[x, j][0] == 
    If[j == 1 && i == 1, 0.7, 
     If[j == i, 
      0.01 Subscript[x, RandomInteger[{1, j - 1}]][t] /. t -> T, 
      Subscript[x, j][t] /. t -> T]]}, {i, n}, {j, i}]

I want solve the differential equation 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).

I'm having trouble building the initial conditions, present in the code above

 Subscript[x, j][0] == 
        If[j == 1 && i == 1, 0.7, 
         If[j == i, 
          0.01 Subscript[x, RandomInteger[{1, j - 1}]][t] /. t -> T, 
          Subscript[x, j][t] /. t -> T]]

Where the first initial condition {x1[0]=0.7}

The next {x1[0]=x1[Last time in j past],x2[0]=[Last time in j past]}

. . .

{{x1[0]=[Last time in j past],x2[0]=x2[Last time in j past],x3[0]==x3[Last time in j past],x4[0]=x4[Last time in j past],x5=0.01xSubscript[RandomInteger[1,j-1]][Last time in j past]}

can anybody help me?

$\endgroup$
6
  • 2
    $\begingroup$ I'm sorry but your question is unclear. What does it mean "ODE system that increases in size according to the rules"? and "This code solve the system of according with interpolation structure {{x1},{x1,x2},{x1,x2,x3},{x1,x2,x3,x4,},{x1,x2,x3,x4,x5}}"? Could you explain what problem you're trying to solve? $\endgroup$ – Chris K Feb 23 '17 at 20:34
  • $\begingroup$ To solve the above system n = 5; T = 50; nu = 0.05; vars = Table[Subscript[x, j][t], {i, n}, {j, i}]; eqns = Table[{Subscript[x, j]'[t] == Subscript[x, j][ t] (1 - Subscript[x, j][t] - nu (Sum[ Subscript[x, k][t] Boole[k != j], {k, i}]) ), Subscript[x, j][0] == If[j == 1 && i == 1, 0.7, If[j == i, 0.01 Subscript[x, RandomInteger[{1, j - 1}]][t] /. t -> T, Subscript[x, j][t] /. t -> T]]}, {i, n}, {j, i}] sol0 = Table[s = NDSolve[eqns[[l]], vars[[l]], {t, 0, T}], {l, n}] I forgot to add the NDSOlve part $\endgroup$ – SAC Feb 23 '17 at 21:31
  • $\begingroup$ @ChrisK The result of interpolation NDSolve provides the results of the form {{x1},{x1,x2},{x1,x2,x3},{x1,x2,x3,x4,},{x1,x2,x3,x4,x5}}. The initial conditions follows the same structure. I need to give the initial condition for each j by taking into account the last value given to j previous. For example, for j=1 Subscript [x, 1] [0] = 0.7. For j = 2, the initial conditions Subscript [x, 1] [0] = Subscript [x, 1] [T], Subscript [x, 2] [0]= Subscript [x, 1] [T] , Where Subscript [x, 1] [T] is the last value of the interpolation obtained at j = 1. $\endgroup$ – SAC Feb 23 '17 at 21:45
  • $\begingroup$ I think verbally describing the problem you want to solve will be clearer than using Mathematica code. $\endgroup$ – Chris K Feb 23 '17 at 21:57
  • $\begingroup$ If you share the mathematical form of the problem then it might help? $\endgroup$ – zhk Feb 24 '17 at 13:15
3
$\begingroup$

This looks like a Lotka-Volterra competition model with weak, symmetric interspecific competition. If I understand correctly, you want to build up a community by introducing new species every T=50 time steps, where the new species has initial density of 0.01 times one of the existing species.

Call me old-fashioned, but I think this is easiest to understand when handled in an iterative Do loop:

nmax = 5; (* max number of species *)
T = 50; (* period *)
nu = 0.05; (* interspecific competition coefficient *)

(* 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]- nu (Sum[Subscript[x, k][t] Boole[k != j], {k, n}]))
, {j, n}];

(* initial ICs *)
ics = {Subscript[x, 1][0] == 0.7};

(* main loop *)
Do[
  (* solve for n species *)
  sol[n] = NDSolve[{eqns, ics}, vars, {t, 0, T}][[1]];
  (* plot dynamics *)
  Print[Plot[Evaluate[Table[Subscript[x, j][t], {j, n}] /. sol[n]], {t, 0, T}, PlotRange -> {0, All}]];

  (* 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, RandomInteger[{1, n}]][T] /. sol[n]]}
   ];
, {n, nmax}]

Mathematica graphics Mathematica graphics Mathematica graphics Mathematica graphics Mathematica graphics

Not sure how interesting this is, since I believe this symmetric LV system has a globally stable equilibrium at $x_i=1/(1-(n-1)nu)$ as long as $0<nu<1$, but you could use this as a basis for more interesting explorations of community assembly.

Edit

Here's a 3D phase portrait of the first three periods for @MMM:

Show[
  ParametricPlot3D[Evaluate[{Subscript[x, 1][t], 0, 0} /. sol[1]], {t, 0, T}],
  ParametricPlot3D[Evaluate[{Subscript[x, 1][t], Subscript[x, 2][t], 0} /. sol[2]], {t, 0, T}],
  ParametricPlot3D[Evaluate[{Subscript[x, 1][t], Subscript[x, 2][t], Subscript[x, 3][t]} /. sol[3]], {t, 0, T}],
  PlotRange -> {{0, 1}, {0, 1}, {0, 1}}
]

Mathematica graphics

$\endgroup$
9
  • $\begingroup$ Thanks for the answer. If possible, please include a phase portrait? $\endgroup$ – zhk Feb 25 '17 at 18:54
  • $\begingroup$ @MMM The dimension of the system keeps increasing, so that'd be tough without ParametricPlot5D! :) :) $\endgroup$ – Chris K Feb 25 '17 at 18:56
  • $\begingroup$ As a test case, restrict it to 3 species or even 2. $\endgroup$ – zhk Feb 25 '17 at 18:59
  • 1
    $\begingroup$ @MMM OK see edit above $\endgroup$ – Chris K Feb 25 '17 at 19:07
  • $\begingroup$ I did. Thankyou $\endgroup$ – zhk Feb 25 '17 at 19:08
1
$\begingroup$

If I understand it correctly, I think the following is what you want:

Solve the differential equation 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).

The following code uses iniCond which generates initial conditions for a given i, diffGen generates equations with the initial conditions, and finally solveDiffGen solves the differential equations.

T = 50;
nu = 0.05;
vars[i_] := Subscript[x, #][t] & /@ Range[i];
iniCond[1] = {Subscript[x, 1][0] == 0.7};
diffGen[i_] := 
  Join[Table[{Subscript[x, j]'[t] == 
       Subscript[x, j][
         t] (1 - Subscript[x, j][t] - 
          nu (Sum[Subscript[x, k][t] Boole[k != j], {k, i}]))}, {j, 
      i}], iniCond[i]] // Flatten;
iniCond[i_] := 
  Module[{prevSolConds, lastCond}, 
   prevSolConds = 
    NDSolve[diffGen[i - 1], 
       vars[i - 1], {t, 0, 
        T}] /. (x_ -> g_) :> {(x /. {t -> 0}) == (g /. {t -> T})} // 
     Flatten;
   lastCond = 
    Subscript[x, i][0] == 
      0.01 Subscript[x, RandomInteger[{1, i - 1}]][t] /. t -> T;
   {prevSolConds, lastCond} // Flatten];
solveDiffGen[i_] := NDSolve[diffGen[i], vars[i], {t, 0, T}];
$\endgroup$
9
  • $\begingroup$ Hi @Anjan this idea is exactly what I'm trying to solve. I tried to use your code suggestion, however the following errors appeared and I could not fix them. Do you know what is happening? $\endgroup$ – SAC Feb 24 '17 at 12:16
  • $\begingroup$ FindRoot::sszero: The step size in the search has become less than the tolerance prescribed by the PrecisionGoal option, but the function value is still greater than the tolerance prescribed by the AccuracyGoal option. >> NDSolve::berr: There are significant errors {-1.,0.} in the boundary value residuals. Returning the best solution found. >> $\endgroup$ – SAC Feb 24 '17 at 12:25
  • $\begingroup$ FindRoot::sszero: The step size in the search has become less than the tolerance prescribed by the PrecisionGoal option, but the function value is still greater than the tolerance prescribed by the AccuracyGoal option. >> NDSolve::berr: There are significant errors {-1.,0.} in the boundary value residuals. Returning the best solution found. >> FindRoot::sszero: The step size in the search has become less than the tolerance prescribed by the PrecisionGoal option, but the function value is still greater than the tolerance prescribed by the AccuracyGoal option. >> $\endgroup$ – SAC Feb 24 '17 at 12:26
  • $\begingroup$ General::stop: Further output of FindRoot::sszero will be suppressed during this calculation. >> NDSolve::berr: There are significant errors {-1.,0.} in the boundary value residuals. Returning the best solution found. >> General::stop: Further output of NDSolve::berr will be suppressed during this calculation. >> $\endgroup$ – SAC Feb 24 '17 at 12:27
  • $\begingroup$ At what value of i is it not working? $\endgroup$ – Anjan Kumar Feb 24 '17 at 12:53

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.